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ABSTRACT 

This  computer  simulation  investigates  the  ways  in  which  an  embedding 
lattice  influences  a  collision  event.   Low  energy  events  were  studied  with 
a  Born-Mayer  potential  function;  high  energy  with  a  Thomas-Fermi-Firsov 
potential.   Conclusions  were:  (1)  the  lattice  will  increase  ranges  for 
particle  energies  above  400  ev;  (2)  the  apparent  mass  concept  is  not  a 
valid  description  of  events  in  a  lattice;  (3)  lattice  effects  will  signi- 
ficantly modify  the  low  energy  portion  of  the  target  atom  energy  distri- 
bution function;  (4)  there  is  no  evidence  that  a  copper  atom  will  "rebound" 
from  a  copper  lattice. 

The  writer  wishes  to  express  his  appreciation  for  the  assistance  and 
encouragement  given  him  by  Associate  Professor  Don  E.  Harrisons  Jr.  of  the 
U.  S.  Naval  Postgraduate  School  in  this  investigation. 
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1.   History. 

Radiation  damage  is  the  overall  term  applied  to  erosion9  disruptions,  or 
rearrangement  of  a  crystal  lattice  by  bombardment  of  ions  or  atoms,,   Sput- 
tering, a  phenomena  sometimes  observed  in  radiation  damage,  applies  to  the 

process  of  lattice  atoms  being  ejected  or  "knocked  out"  of  the  crystal.   The 

1 
sputtering  process  was  first  observed  in  1852  by  Grove  in  the  cathodes  of 

electrical  discharge  tubes. 

Since  Grove's  observation,  experimental  and  theoretical  physicists  have 

been  attempting  to  develop  a  sputtering  theory  that  will  predict  measurable 

experimental  quantities  as  a  function  of  mass  ratio,,  position  in  the  periodic 

table,  energy  of  bombarding  particle  and  angle  of  incidence.   To  date9  no 

such  complete  theory  is  available. 

A  tremendous  volume  of  material  has  been  written  on  the  subject  of  spuit- 

2 
tering  since  its  discovery.   In  1955,  Irvine  listed  306  papers  concerned 

with  the  sputtering  process.   These  early  papers  disagree  violently  because 

important  experimental  conditions  were  not  obtained  or  simply  ignored; 

consequently  many  of  the  conclusions  reached  are  contradictory. 

During  this  period  two  important  theories  were  proposed.   Kingdon 

3 
Langmuir  ,  in  1923,  developed  a  momentum  transfer  theory  which  agreed 

well  with  experimental  data  they  had  obtained  for  thoriated  tungsten  fila* 

ments.  The  "evaporation  theory"  was  proposed  in  1928  by  von  Hippel  and 

4  5 

Blechschmidt  and  improved  by  Townes   in  1944.   It  assumes  that  the  inc 


W.  R.  Grove,  Phil.  Trans.  Roy.  Soc.  London,  142,  87  (1842). 

2 
M.  M.  Irvine,  Dissertation,  Lehigh  University  (1955). 

3 
K.  H.  Kingdon  and  I.  Langmuir,  Phys.  Rev.  22,  148  (1923). 

4 
A.  von  Hippel  and  E.  Blechschmidt,  Ann.  Physik  86,  1006  (1928) 

C.  H.  Townes,  Phys.  Rev.  65,  319  (1944). 


particles  heat  a  small  volume  of  the  crystal  to  a  very  high  temperature 

6 
and  subsequently  surface  atoms  evaporate.   More  recently  (1953) 9  Harrison  9 

following  a  suggestion  made  by  Keywell  s  presented  a  statistical  theory  of 
sputtering  which  could  be  placed  somewhere  between  the  two  models  mentioned 
above.  We  will  not  present  details  of  these  theories.   We  must  mention  how- 
ever, that  they  do  not  provide  an  adequate  description  of  the  sputtering 
process. 

The  trend  in  sputtering  theories  has  recently  been  toward  the  considers 

8 
tion  of  the  individual  collision  processes  within  the  lattice.  Henschke 

i 
has  used  this  approach  in  his  theory  of  sputtering  and  Harrison  and  Magnuson' 

have  also  applied  it,  with  modifications,  to  the  theoretical  study  of  sput- 
tering threshold  energies. 

The  theoretical  work  performed  by  Henschke  uses  the  momentum  transfer 

3 
concept  initiated  by  Kingdon  and  Langmuir  to  account  for  all  the  experi- 

mental  phenomena  observed  at  low  incident  particle  energy.   His  basic  treat- 

ment  of  collisions  between  atoms  in  a  lattice  is  based  on  the  assumption 


...that  the  collisions  described  can  be  treated  with  the  general 
principles  of  classical  mechanics,  using  impulsive  forces,  in  a 
manner  similar  to  the  well-known  collisions  with  restitution.  " 


However,  he  did  not  in  all  cases  use  the  masses  of  the  individual  particles 
involved  in  the  equations  pertaining  to  the  two  body  collisions.   Instead 
he  postulated  an  "effective"  mass.   If  the  collision  is  between  a  surface 


D.  E.  Harrison,  Jr.,  Phys.  Rev.  102,  1473  (1956). 

7F.  Keywell,  Phys.  Rev.  87,  160  (1952). 

8E.  B.  Henschke,  Phys.  Rev.  106,  737  (1957). 
9 
D.  E.  Harrison  and  G.  D.  Magnuson,  Phys.  Rev.  122,  1421  (1961) 

See  reference  8,  p.  738. 


atom  and  a  moving  ion  or  atom  and  is  such  that  the  surface  atom  is  st 
on  its  "inside"  hemisphere,  according  to  Henschke,  an  effective  mass  is 
required  and  the  two  particles  involved  can  be  considered  to  have  their 
masses.   Should  the  collision  of  an  atom  or  ion  be  directed  inward  from  the 
target  lattice  surface,  Henschke  states, 


The  bulk  of  the  target  is  behind  the  struck  atom  and  produces  a 
very  large  'effective'  mass  compared  to  the  mass  of  the  ion  or 
to  the  mass  of  the  target  atom.H 


The  reason  Henschke  did  not  assume  pure  elastic  collisons  was  stated  as 


Energy  losses  are  due  to  the  fact  that  the  target  atom  is  coupled 
rather  strongly  with  the  atoms  of  the  lattice.   Before  the  moment 
of  highest  compression  is  reached,  the  ion  and  the  struck  atom 
exchange  energy  with  the  neighboring  atoms  of  the  lattice.   Debye 
waves  are  thus  excited  and  dissipated  irreversibly  into  the 
tice.12 


He  goes  on  to  state, 


...the  final  step  in  each  sputtering  process  at  any  angle  of  inci- 
dence of  the  ion  can  be  generally  described  as  a  collision  of  the 
ion  with  an  upper  surface  atom,  in  which  this  atom  is  hit  on  its 
inside  hemisphere  so  as  to  obtain  an  impulse  with  a  component  in 
the  direction  of  the  outward  normal  to  the  surface.   If  the  energy 
transferred  in  this  direction  to  this  target  atom  by  the  impart 
of  the  ion  is  equal  to  or  greater  than  the  heat  of  vaporization^ 
with  which  the  atom  is  assumed  to  be  bound  to  the  crystal  lattice 
plane,  then  this  atom  is  ejected  in  the  collision.  J 


Early  study  in  radiation  damage  was  precipited  by  the  advent  of 
nuclear  reactor.   Recently  however,  with  the  imminent  possibility  of  thermo- 
nuclear power  production  and  ion  propulsion  engines  for  space  vehicless  othes 
practical  applications  of  the  processes  involved  in  radiation  damage  have 


11E.  B.  Henschke,  Phys.  Rev.  121,  1290  (1961). 


12 

See  reference  8,  p.  738. 

13 

See  reference  8,  p.  737. 


become  extremely  important. 

14 
Robinson  and  his  co-workers   have  made  theoretical  studies  of  the 

ranges  in  solids  of  atoms  having  energies  from  1  to  100  kev  using  digital 
computer  techniques  on  the  basic  assumption  that  the  moving  atom  loses  its 
energy  through  repeated  binary  elastic  collisions  with  atoms  of  the  solid. 
The  masses  used  in  these  calculations  are  the  true  masses  of  the  interact- 
ing particles,  but  they  note  that  their  assumption  is  certainly  not  valid 

15 
below  100  ev.   Gibson  et.  al.   have  also  taken  advantage  of  the  speed  avail- 


able with  modern  digital  computers  in  the  study  of  radiation 
do  not  assume  binary  collisions  but  instead  use  iteration  techniques  eimp] 
ing  Newton's  equation  of  motion  to  solve  the  complex  many  body  problem.  The 
initial  success  of  these  programs  have  encouraged  their  originators  to  ex- 
plore further  the  possibilities  of  these  techniques.   Although  the  computer 
programs  are  designed  to  study  radiation  damage,  the  basic  principles  involv- 
ed are  also  important  in  sputtering. 

The  work  just  discussed  is  very  sensitive  to  the  mathematical  form  of 
the  interatomic  potential.   The  Born-Mayer,  Bohr  (screened  Coulomb) 9  and 
the  Thomas-Fermi-Firsov  are  some  of  the  potentials  used  in  these  calcula- 
tions.  The  ability  of  any  of  the  potential  functions  to  describe  a  physical 
situation  is  highly  dependent  upon  the  energy  range  under  consideration  and 
the  specific  atoms  or  ions  involved.   No  single  potential  has  yet  been  de- 
vised that  satisfactorily  represents  the  Intel  action  ander  all  circumstances, 


14 

0.  S.  Oen,  D.  K.  Holmes,  and  M.  T.  Robinson,  Jour  Appl.  Phys.  34  s 

302  (1963). 

J.  B.  Gibson,  A.  N.  Goland,  M.  Milgram,  and  G.  H.  Vineyard,  Phys, 


Rev.  120,  1229  (1960) 


An  obvious  unresolved  question  is  inherent  in  these  theories „  The 
many  body  approach  and  the  two  body  assumption  with  a  minimum  energy  limita- 
tion, appear  to  be  reasonable  approaches  to  the  problem.,   Both  methods  des- 
cribe possible  events  in  a  lattice,  but  are  they  compatable?  The  general 
use  of  the  two  body  approach  is  desirable  because  of  its  simplicity,,   The 

problem  to  be  studied  is:  when  can  the  two  body  collision  be  &ssumed9  if 

a 

at  all,  and  should  an  "apparent1'  mass,  as  proposed  hy 


2.   Objective. 

We  intend  to  investigate  theoretically  the  energy  transfer  process 

between  an  incoming  copper  atom  and  a  copper  lattice.  The  many  body  ap- 

15 
proach,  using  computer  techniques  similar  to  Gibson  et.  al.s  will  be 

utilized  and  then  compared  with  the  simple  two  body  solution.   We  hope 

the  results  will  provide  some  quantitative  answers  to  the  questions  posed 

in  the  preceding  section. 


3.   Procedure. 

The  representative  crystal  lattice  consists  of  a  cube  of  63  atoms 
arranged  in  a  face  centered  cubic  structure 0   This  arrangement  gives  a 
total  of  63  atoms o  These  atoms  are  free  to  move,  if  disturbed,,   In  an 
effort  to  simulate  a  larger  lattice,  the  movable  63  atoms  are  surrounded 
on  all  sides  by  stationary  immovable  atoms  positioned  as  a  continuation  of 
the  crystal.  This  procedure  also  gives  the  bombarded  face  some  character- 
istics similar  to  a  binding  energy. 

In  order  to  simplify  the  calculations,  we  defined  a  quantity  called 
the  "lattice  unit".   One  lattice  unit  is  equal  to  half  the  length  of  a 
cube  side.   For  copper,  the  length  of  a  cube  side  is  3.614A,  therefore 
one  lattice  unit  is  1.807A.  A  three  dimensional  Cartesian  coordinate 
system  is  used  for  reference  with  the  origin  placed  with  a  (100)  face  of 
the  movable  63  atom  "core"  parallel  to  the  x-z  plane  at  y=0.  Any  penetra- 
tion into  the  lattice  core  from  this  plane  is  in  the  +  y  direction.  The 
immovable  atoms  surrounding  the  core  are  in  the  x-z  plane  at  y=-l  and  y^5, 
the  x-y  plane  at  z=0  and  zs=6,  and  the  y-z  plane  at  x^O  and  x^6.   For  convem- 
ience  all  atoms  are  numbered.   Number  one  is  the  "bullet'0,  numbers  2-64  are 
the  "core"  and  65-172  are  the  immovable  atoms. (see  Fig.  1) 

The  program  depends  upon  a  potential  function  subroutine;  so  potential 
modification  is  accomplished  without  any  changes  in  the  main  program.  The 
Born-Mayer  and  the  Thomas-Fermi-Firsov  have  been  used  in  the  present  work 
(see  Appendix  VI  for  exact  forms  and  constants).  All  calculations  are 
with  the  potential  in  "eroded"  form.  To  "erode"  a  potential,  the  value  of 
the  potential  at  the  nearest  neighbor  separation  is  subtracted  from  the 
value  of  the  potential  at  smaller  distances.   The  eroded  form  makes  the 


potential  zero  at  distances  greater  than  nearest  neighbors. 

The  interaction  between  atoms  is  represented  by  an  eroded  repulsive 
potential.  Also,  the  entire  lattice  is  initially  at  absolute  zero  since 
no  vibrational  energy  is  simulated.  These  two  approximations  imply 
the  lattice  has  no  potential  or  kinetic  energy  before  interaction  with 
bullet  and  all  the  energy  in  the  lattice  at  any  time  thereafter  is  deri 
from  the  bombarding  atom. 

We  realize  that  each  atom  in  the  lattice  is  contained  in  some  type  of 
"potential  well",  but  we  are  measuring  the  energy  from  the  ground  state 
level  of  this  "well"  rather  than  from  the  true  zero  of  potential  at  infinite 

separations.  Thus,  before  interactions  with  the  bullet,  our  lattice  is  not 

15 
held  together  by  the  surrounding  stationary  atoms  as  in  the  Gibson  et.  al. 

model  because  the  erosion  also  has  removed  all  forces.   Forces  appear  only 

when  atoms  move  from  their  equilibrium  positions. 

The  use  of  an  eroded  potential  does  not  provide  as  good  a  model  as  the 
one  used  by  Gibson  et.  al.,  but  it  reduces  the  computer  time  required  for 
calculations. 

The  "binding  energy"  created  by  the  layer  of  fixed  atoms  covering  the 
front  face  is  probably  much  too  large,  but  it  does  approximately  simulate 

the  behavior  of  the  next  layer  of  atoms. 

1  ^ 
The  calculations  which  move  atoms  are  similar  to  those  of  Gibson  et.  al. 

Newton's  equation  of  motion  can  be  rearranged  to  give  the  change  in  velocity 
of  a  body  acted  on  by  an  unbalanced  force.   (     F  '^ ■'/ pTl  ^  ^— ^    ) 
The  change  in  velocity  can  be  related  to  a  change  in  position  if  an  average 
velocity  is  assumed  (  Xjjew  "[V ^"^"j  ^   "*  ^N^q).  Tne  unbalanced  force 
used  in  the  first  equation  is  an  average  force  calculated  by  a  double  itera- 
tion procedure  as  follows;   (1)  assume  an  atom  at  position  1  with  velocity  1 

8 


(2)  calculate  the  total  force  on  the  atom  as  a  result  of  all  the  other  atoms 
in  the  lattice  (this  means  normally  only  about  8-10  nearest  atoms  because 
the  potential  is  eroded)   (3)  call  this  calculated  forces  force  ls  and  use  the 
equation  of  motion  to  move  the  atom  to  a  temporary  posit ion8  position  2„ 

(4)  now  repeat  the  force  calculations  for  position  2,  call  this  force  20 

(5)  go  back  to  position  1,  and  use  the  average  of  force  1  amid  force  2  to 
move  the  atom  to  a  new  position,  position  3.   Procedures  1  through  5  consti- 
tute one  "time  step!'.   Forces  are  eroded  in  the  same  manner  as  potentials 
and  are  calculated  by  a  subroutine  based  on  the  partial  derivative  of  the 
potential  function  with  respect  to  distance. 

The  basic  value  of  AT  is  dependent  upon  the  original  energy  of  the 
bombarding  atom.   A  T  is  the  time  in  seconds  required  for  the  incoming  atom 
to  traverse  one  lattice  unit.  The  atom  will  lose  energy  as  it  interacts 
with  the  lattice  but  AT  remains  unchanged.   To  add  flexibility  to  the 
program,  the  basic  value  of  AT  can  be  multiplied  by  any  desired  factor 
(time  step  multiplier).   A  T  is  established  at  the  beginning  of  the  program 
and  is  constant  for  the  duration  of  the  calculations. 

The  bombarding  copper  atom  is  originally  located  in  the  plane y  ■  -v2 
and  given  a  velocity  in  the  +  y  direction  only.  This  requires  the  bombarding 
particle  to  approach  perpendicular  to  the  (100)  face  of  the  movable  lattice 
core.   The  immovable  atoms  in  the  plane  (y  =  -1)  that  cover  the  bombarded 
face  do  not  interact  with  the  bullet. 

An  impact  area  was  chosen  on  the  (100)  face  so  that  all  points  in  this 
area  would  be  representative  of  any  point  in  a  (100)  plane.  An  impact  point 
is  defined  as  the  location  on  the  (100)  face  toward  which  the  bombarding  atom 
is  directed.   It  may  or  may  not  actually  pass  through  this  point.   By  moving 
in  1/10  lattice  unit  increments,  36  impact  points  in  the  impact  "triangle" 


were  assigned.   (see  Fig.  2  and  3). 

For  each  case  the  two  body  problem  is  solved  after  the  many  body 
problem  with  the  same  double  iteration  scheme.   Lattice  atom  number  eight 
was  assigned  as  the  "target"  in  the  many  body  lattice  problem  so  that  its 
motion  could  be  compared  with  the  simple  two  body  problem.  The  geometrical 
relationship  and  physical  constants  of  the  two  body  problem  correspond  to 
those  of  the  incoming  atom  and  "target"  in  the  crystal. 

The  program  was  written  in  FORTRAN  language  for  use  on  a  CDC  1604 
computer.  The  input  consists  of  the  incoming  particle  mass  and  energy  and 
the  time  step  multiplier.   From  these  inputs  the  bullet  is  assigned  a  loca- 
tion and  an  appropriate  velocity  component  in  the  +  y  direction.   Also  cal- 
culated is  the  actual  value  of  AT. 

The  next  section  of  the  program  assigns  coordinates  to  all  atoms  in 
the  crystal,  movable  and  immovable,  gives  each  atom  a  number  and  sets  their 
velocities  equal  to  zero.   If  the  distance  between  any  two  atoms  is  less 
v2  lattice  units,  the  force  between  these  two  atoms  is  calculated  and  stored 
in  memory.   This  is  done  for  each  atom  and  the  summation  of  all  the  forces 
calculated  gives  the  resultant  force  on  each  particle. 

The  change  in  the  velocities  of  all  movable  atoms  as  a  result  of  un- 
balanced forces  is  calculated;  then  the  original  velocities  of  all  particles 
and  the  calculated  change  in  velocity  are  used  to  move  atoms  to  temporary 
locations.   The  calculations  to  obtain  resultant  forces  are  repeated  for  the 
temporary  locations,  old  and  new  forces  are  averaged  and  used  to  move  all 
atoms  to  new  locations.   Forces  and  velocities  are  calculated  as  vector 
components  using  the  coordinate  system  established  for  the  lattice. 

After  each  movement  of  all  the  atoms,  the  kinetic  energy  of  the  tar- 
get is  calculated  and  compared  to  its  former  value.  The  procedure  is  repeated 
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until  the  kinetic  energy  of  the  target  reaches  a  maximum.  The  potential 
and  kinetic  energy  of  each  movable  atom  is  then  calculated.  The  total 
energy  of  the  entire  lattice  including  bullet  is  also  found  and  compared  to 
the  original  "bullet"  energy.  This  comparison  gives  a  reasonable  measure  of 
the  accuracy  achieved.  The  two  body  problem  is  solved  in  the  next  section 
using  the  same  iteration  technique.  Appendices  1-7  contain  a  detailed  ex= 
planation  of  the  program. 

For  each  "run",  the  computer  will  solve  the  many  body  and  two  body 
problem  36  times,  corresponding  to  the  36  different  impact  points.  The  form 
of  the  printed  output  for  each  impact  point  is  shown  in  Fig.  4  and  5.   The 
"triangular"  output  as  seen  in  Fig.  6  indicates  how  two  chosen  parameters 
vary  over  the  impact  area. 

The  computer  running  time  for  each  of  the  36  impact  points  was  2-4 
minutes,  the  time  varying  with  the  value  chosen  for  AT.  We  origin 
assumed  that  the  accuracy  of  the  results  would  increase  for  smaller 
smaller  values  of  AT.  We  found  that  the  value  of  AT  which  will  give  the 
minimum  error  in  total  energy  appears  to  be  a  complicated  function  of  the  time 
step  multiplier,  the  incoming  bullet  energy,  and  the  impact  point.   By  trial 
and  error  methods,  we  found  that  below  100  ev  for  the  Born-Mayer  potentials, 
a  time  step  multiplier  of  around  0.09  would  produce  reasonable  errors  at  all 
impact  points. 
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4.  Potentials 

The  Born-Mayer  potential  is  an  exponential  function  of  the  internuclear 
separation  with  constants  which  may  be  obtained  from  experimentally  measur- 
able elastic  moduli.   The  constants  used  in  the  program  are  the  same  as 
those  of  Gibson  et.  al.   ,  with  special  emphasis  on  their  potential  number 
two.   Near  equilibrium  separations,  the  Born-Mayer  potential  is  thought  to 
be  an  adequate  approximation.  We  expect  the  potential  to  fail  for  energies 
above  a  few  hundred  electron  volts. 

The  Thomas-Fermi-Firsov  (TFF)  potential  is  the  result  of  theoretical 

17 
work  by  Firsov   based  on  the  Thomas-Fermi  model  of  the  atom.   It  is  a 

screened  Coulomb  potential  with  a  more  complicated  screening  function  than 

that  used  in  the  Bohr  potential.  The  TFF  is  satisfactory  for  intermediate 

separations.  This  places  it  between  the  Bohr  potential  used  for  small  separa- 

18 
tions  and  the  Born-Mayer  potential.   Abrahamson  and  Hatcher   state  that  the 

Thomas-Fermi  approximation  becomes  unreliable  when  the  internuclsar  distance 

exceeds  f\J  1A.   The  form  of  the  TFF  potential  and  the  appropriate  constants 

used  in  the  program  are  the  same  as  those  used  in  computer  programs  at  the 

19 

Oak  Ridge  National  Laboratory 

The  mathematical  forms  of  the  potential  and  force  functions  as  they  ap- 
pear in  the  program  are  derived  in  Appendix  VI, 


See  reference  15,  p.  1233. 

170.  B.  Firsov,  J.  Exptl.  Theoret.  Phys.  (D.S.S.R.)  328  1464  (1957), 

translation:  Soviet  Phys.  -JETP  5,  1192  (1957). 

1  ft 

A.  A.  Abrahamson  and  R.  D.  Hatcher,  Phys.  Rev.  121,  159  (1961). 

19 

R.  T.  Robinson  (private  communication).  The  author  would  like  to 

thank  Dr.  Robinson  for  the  use  of  his  results  prior  to  their  publication. 
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5.   Results 

A.  General  description  of  events  In  the  lattice. 

The  number  of  atoms  disturbed  in  the  lattice  is  a  function  of  two 
variables,  the  number  of  time  steps  the  computer  performs  and  the  impact 
point  of  the  bullet.   In  general,  the  lattice  is  more  disrupted  by  impact 
points  in  the  area  near  x  =  3.0,   z  ■  4.0,  and  the  lattice  is  least  af- 
fected by  head  on  collisions  with  the  target.  No  more  than  one  half  of 
atoms  were  ever  disturbed  before  the  target  reached  a  maximum  in  kinetic 
energy. 

The  impact  area  is  not  symmetrically  located  with  respect  to  the  sur- 
rounding immovable  atoms  but  for  appropriate  points  in  the  impact  area8 
excellent  symmetry  of  displacement  was  observed  for  the  movable  atoms.  This 
indicates  that  the  size  of  the  model  is  adequate  for  present  purposes. 

Because  of  the  limitation  placed  on  the  number  of  time  steps ,  only  a 
few  atoms  have  any  significant  kinetic  energy  when  the  interaction  is  stop- 
ped. These  are  normally  the  target,  bullet,  and  the  atoms  directly  behind 
the  target  and  bullet  in  the  lattice. 

Although  our  primary  purpose  was  to  compare  lattice  interactions  with 
the  two  body  interaction,  we  altered  the  program  so  that  it  would  run  a 
predetermined  number  of  time  steps  in  order  to  observe  "chains"  and  "chan- 
nels".  Impact  points  at  x  =  3.0,  z  =  4.0  and  x  ■  3.0,  z  ■  3.0  correspond 
to  the  beginning  of  (100)chains.  These  are  not  "close  packed"  chains  but 
they  transfer  energy  reasonably  well,  especially  above  100  ev.   Impact  point 

x  ■  2.5,  z  m   3.5  is  the  beginning  of  a  (100)  channel.  The  (100)  channel  is 

.* — 
not  as  "wide  open"  as  a  (110)  channel  and  a  100  ev  bullet  does  not  travel 

more  than  f*J  3A  into  the  lattice. 
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B.  Comparison  of  lattice  interactions  with  the  two  body  interaction. 

A  strict  comparison  of  interactions  in  the  lattice  to  the  two  body 
problem  is  not  possible,  but  certain  limited  correlations  can  be  made, 
difficulty  arises  because  the  end  of  the  interaction  in  the  lattice  cannot 
be  defined  explicitly.  There  are  several  criteria  available  that  may  be 
used  in  an  attempt  to  compare  the  interactions.   If  the  two  body  interaction 
was  allowed  to  continue  until  the  potential  between  atoms  was  zero8  for 
comparison  purposes  the  lattice  interaction  could  procede  until;  1)  the 
potential  between  target  and  bullet  is  equal  to  zero  2)  the  total  energy 
(kinetic  plus  potential)  of  the  target  reaches  a  maximum  or  3)  the  kinetic 
energy  of  the  target  reaches  a  maximum.   All  of  the  criteria  apply  to  the 
same  physical  situation  in  the  two  body  problem  but  are  not  equivalent  in 
the  lattice  interaction. 

We  chose  the  maximum  kinetic  energy  of  the  target  as  the  comparison 
stopping  point  because  this  represents  the  time  when  the  target  and  bullet 
have  almost  ceased  to  interact  and  are  just  beginning  to  interact  with  other 
lattice  atoms. 

A  fourth  criteria  was  suggested  after  the  present  results  had  been 
obtained.   When  the  kinetic  energy  of  the  target  in  the  lattice  reaches  a 
maximum,  find  the  distance  between  the  target  and  bullet  and  allow  the  two 
body  interaction  to  proceed  until  the  distance  between  atoms  is  equal  to  the 
separation  in  the  lattice.  This  criteria  is  possibly  better  than  the  one 
chosen  but  at  the  present  time  it  has  not  been  utilized. 

The  geometrical  relationships  of  the  target,  bullet,  and  lattice  intro- 
duce complications  into  the  comparison  scheme.   The  motion  of  the  target  in 
the  lattice  is  not  restricted  isotropically   (see  Fig.  7).  If  the  inter- 
action with  the  bullet  moves  the  target  toward  atom  number  six  (shaded  area)s 
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the  movement  of  the  target  is  restricted  more  than  if  the  target  motion  is 
initially  toward  the  cross-hatched  area.   Restriction  of  the  target  motion 
lengthens  the  interaction  time  between  bullet  and  target,  therefore  the 
energy  transferred  to  the  target  is  greater  than  in  the  two  body  problem. 
This  effect  is  noticable  only  at  impact  parameters  of  0.5  lattice  units  or 
greater.   At  impact  distances  smaller  than  this  the  interaction  occurs  so 
rapidly  that  the  target  does  not  move  a  significant  amount  and  the  effect  is 
reduced.   The  kinetic  energy  transferred  at  distances  smaller  than  0.5  lattice 
units  is  less  than  that  transferred  in  the  two  body  case  because  there  is 
absorption  of  energy  by  the  remainder  of  the  lattice.  The  geometry  of  the 
bombarded  face  is  not  the  only  factor  that  effects  energy  transfer.  The 
atoms  behind  the  target  in  the  lattice  also  play  an  important  role  but  their 
effect  is  not  as  immediately  obvious  as  those  just  discussed. 

Figure  8  illustrates  some  of  the  points  just  mentioned.   Along  the  line 
A  B  (x  =  3.0,  z  =  3.0-4.0)  the  target  behavior  in  the  lattice  is  similar  to 
the  two  body  interaction.  These  impact  points  cause  the  target  to  move  in 
the  negative  z  direction  toward  a  "hole"  in  the  lattice  face.   At  point  D  on 
the  line  A  C  ,  the  effect  of  restricted  target  motion  becomes  obvious,  and 
also  in  the  area  enclosed  by  the  dashed  line.   The  geometrical  effects  of 
the  lattice  are  fairly  predictable  and  no  unusual  phenomena  are  observed. 

For  a  static  system  of  two  equal  mass  atoms  with  a  conservative  repul- 
sive force  between  them,  the  total  potential  energy  of  the  system  can  be 
halved  and  the  result  assigned  as  the  energy  of  the  atom.   If  the  two  atoms 
are  allowed  to  move,  the  total  energy  (potential  plus  kinetic)  of  each  will 
be  constant  and  equal  to  the  original  potential  energy.   In  a  system  composed 
of  three  or  more  atoms,  this  is  generally  not  true.   For  any  system  where  the 
atoms  are  not  static,  the  process  of  assigning  an  energy  to  any  one  atom  is 
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no  longer  possible  (there  is  one  exception).   This  is  the  situation  in  the 
lattice. 

The  potential  energy  of  an  atom  in  the  lattice  is  defined  in  the  program 
as  half  the  potential  associated  with  its  position.  This  definition  is  merely 
a  convenience. 

The  area  of  the  impact  triangle  from  point  A  to  the  arc  E-F  (see  Fig.  8) 
is  least  effected  by  the  geometry  of  the  lattice  and  closer  inspection  of 
these  interactions  is  justified.   The  kinetic  energy  transferred  to  the 
n-body  target  is  less  than  that  of  the  two  body  target.   In  percentages8 
the  kinetic  energy  transferred  becomes  less  as  the  bullet  energy  is  decreased. 
(see  Figs  8  through  11).   In  each  case  the  target  also  acquires  some  poten- 
tial energy  (see  Figs.  12  through  15).   As  stated  earlier,  it  is  difficult  to 
assign  a  specific  potential  energy  to  the  target.   We  do  know  the  limits 
involved,  i.e.  the  target  could  eventually  receive  none  of  it,  or  twice  the 
indicated  value.   Either  limiting  situation  is  unlikely. 

If  we  make  the  assumption  that  one-half  of  the  potential  "belongs" 
to  the  target,  then  the  total  energy  of  the  target  is  approximately  equal  to 
the  energy  transferred  in  the  two  body  problem  (within  27»)  for  energies  great- 
er than  50  ev. 

Bullet  behavior  for  the  area  under  discussion  is  also  very  similar  to 
the  two  body  problem  (see  Figs.  16  through  19).   For  the  remainder  of  the 
impact  triangle,  the  bullet  is  strongly  effected  by  other  atoms  in  the  lat- 
tice. 

Three  runs  were  made  with  the  TFF  potential  at  energies  of  5,  10s  and 
30  kev  (Figs.  20  through  25).   Since  the  TFF  potential  is  unreliable  beyond 
1A,  the  section  ABC  in  Figs.  20  through  25  is  the  only  part  of  the  impact 
area  that  can  be  considered.   This  section  includes  all  impact  points  that 
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result  in  a  closest  point  of  approach  of  1A  or  less.   We  noted  the  same 
general  agreement  with  the  two  body  interaction  as  found  for  the  Born-Mayer 
potential.  The  only  apparent  difference  is  a  shift  in  the  energy  scale. 

For  a  head-on  collision  in  the  lattice  there  is  very  little  geometrical 
effect,  and  an  "effective"  mass  can  be  described.  At  25  ev,  the  mass  of  the 
target  is  apparently  about  2.5  times  the  bullet  mass,  but  the  energy  trans- 
fer  is  still  approximately  80%.   For  a  500  ev  bullet,  the  effective  mass 
increase  is  negligible.  We  used  the  program,  as  modified  to  observe  chains 

and  channels,  and  could  find  no  evidence  of  bullet  recoil  which  could  be 

8 
attributed  to  an  effectively  heavy  target  as  proposed  by  Henschke  .  Although 

an  effective  mass  concept  for  the  head-on  collision  is  possible,  it  does  not 
properly  describe  the  subsequent  motion  of  either  target  or  bullet.   If  the 
bullet  motion  results  in  a  glancing  hit  with  the  target,  no  single  effective 
mass  can  be  assigned  to  the  target  because  it  is  a  function  of  the  impact 
parameter.   An  average  effective  mass  is  not  applicable  because  the  geometri- 
cal effects  of  the  lattice  on  the  direction  and  energy  of  the  recoil  atom  are 
far  more  significant  than  the  mass  of  the  target  in  a  glancing  hit. 

These  results  apply  to  a  collision  anywhere  in  the  lattice  and  are  not 

limited  to  the  interaction  of  surface  atoms  with  incoming  particles.   The 

g 
results  indicate  that  Henschke' s  apparent  mass  concept  and  the  rebound 

phenomena  associated  with  it  are  not  a  good  description  of  collision  events 

in  a  lattice. 
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6.   Conclusions 

A.  The  Bullet. 

For  energies  above  40  ev  in  the  impact  area  AEF  (see  Figs  16  and 
the  n-body  bullet's  kinetic  energy  after  interaction  is  essentually  the 
same  as  the  two  body  bullet.   Between  40  and  400  ev  for  impact  points  out- 
side the  area  AEF,  the  n-body  bullet  has  considerably  less  energy  after  inter- 
action than  the  two  body  bullet.   Above  400  ev,  the  n-body  and  two  body 
energy  transfers  agree  within  3%  for  all  points  in  the  impact  area. 

The  angular  behavior  of  a  bullet  in  the  lattice  is  much  more  complicated. 
Even  at  energies  above  400  ev,  the  scattering  angle  in  the  n-body  problem  is 
affected  by  the  geometry  of  the  lattice.   In  general,  the  scattering  angle  in 
the  lattice  is  smaller  than  that  of  the  two  body  problem.   It  is  possible  that 
the  lattice  is  attempting  to  focus  the  energy  into  preferred  directions.   Our 
bullets  in  the  (100)  direction  appear  to  focus  in  the  (100)  directions  especial 
ly  at  higher  energies. 

The  cummulative  effects  of  energy  transfer  and  lattice  geometry  will 
affect  the  range  of  energetic  atoms.   The  reduced  scattering  angle  in  the 
n-body  model  for  atoms  with  energies  above  400  ev  should  lead  to  ranges  that 
exceed  those  found  by  two  body  approximation  methods.   Below  400  evs  the 
bullet  loses  more  energy  than  in  the  two  body  case;  so  the  effects  may  cancel, 
or  perhaps  ranges  calculated  by  two  body  methods  may  be  too  large. 

B.  The  Target. 

In  the  area  AEF  (see  Figs.  8,  9,  12,13),  the  n-body  energy  transferred 
to  the  target  is  the  same  (within  27»)  as  that  in  the  two  body  problem  for 
energies  above  50  ev.   For  the  remainder  of  the  impact  area,  the  n-body 
target  receives  more  energy  than  the  two  body  target  for  all  energies  above 
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25  ev.   Isolated  points  at  various  energies  may  not  conform  to  this  rule  but 
the  effect  is  always  present  for  impact  points  near  the  apex  of  the  impact 
area  (x=2.5,  z=3.5).   Near  the  apex,  if  the  bullet  energy  is  large  (300  ev 
or  greater),  the  energy  transferred  to  the  target  is  a  very  small  fraction 
of  the  total  for  both  the  n-body  and  two  body  problems.   However,  the  energy 
transferred  to  the  target  in  the  two  cases  can  differ  by  as  much  as  75%  of  the 
transferred  energy. 

As  a  consequence  of  this  phenomena,  the  low  energy  portion  of  the  n- 
body  target  atom  energy  distribution  (assuming  more  than  one  collision  has 
occurred)  will  be  much  higher  than  that  expected  with  a  two  body  collision 
assumption  for  lattice  interactions.   Precise  measurements  of  energy  trans- 
fer for  large  impact  parameters  do  not  appear  to  be  warranted  in  view  of  the 
results  obtained. 

Unfortunately,  the  face  centered  cubic  lattice  structure  was  incorporat- 

20 
ed  in  the  program  before  we  learned  of  Veksler's   experimental  work  with 

molybdenum  targets  (bcc),  but  certain  qualitative  comparisons  are  still 
possible.   Veksler  has  interpreted  the  lattice  behavior  in  terms  of  the  ef- 
fective mass  concept  but  his  general  conclusions  are  consistent  with  our 
work.   Our  model  gives  good  evidence  to  support  Veksler's  position  that  the 
pair  collision  model  using  elastic  spheres  is  not  acceptable. 

Two  body  approximations  inherently  imply  that  certain  information  is 
available  from  preliminary  n-body  calculations.  This  work  attempts  to  answer 
certain  questions  about  interactions  and  cross  sections.   As  anticipated,,  the 
results  are  not  definitive  but  they  do  indicate  sensitive  areas  which  require 
further  examination. 

20 

V.  I.  Veksler,  J.  Exptl.  Theoret.  Phys.  (U.S.S.R.)  42,  325  (1962). 

translation:  Soviet  Physics  -  JETP  15,  222  (1962). 
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Fig.  1.    THE  MODEL 
Large  circles  indicate  movable  atoms; 
small  circles  indicate  fixed  atoms. 
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Fig.  2.  FACE  (100)  of  LATTICE  CORE. 
Numbers  indicate  atoms;  shaded  portion  is  impact  area 
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Fig.  3.  IMPACT  AREA 
Numbers  indicate  impact  points. 
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Fig.  7.    GEOMETRICAL  EFFECTS  OF  THE  LATTICE  FACE. 
Target  movement  is  more  restricted  if  initial  motion  is  in  shaded  area 
rather  than  in  crosshatched  area. 
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APPENDIX  I 
GENERAL  DISCUSSION  OF  THE  PROGRAM 

The  masses  of  the  lattice  atoms  are  equal  to  that  of  copper  and  are 
an  integral  part  of  the  program.   The  bullet  mass  however,  can  be  varied  by 
changing  the  input  data.  The  potentials  and  forces  are  calculated  by  func- 
tion subroutines  at  the  end  of  the  program.   The  main  program  is  not  dis- 
rupted when  a  different  potential  function  is  substituted. 

A  coordinate  system  was  established  with  an  atom  located  at  the  origin. 
For  a  face  centered  cubic  lattice  in  this  coordinate  system,  the  sum  of  the 
coordinates  of  an  atom  is  always  an  even  number.   Fixed  point  numbers  were 
used  for  lattice  units,  and  the  volume  which  contains  all  the  atoms  (movable 
and  fixed)  was  scanned  in  the  x,  y,  and  z  directions.   If  the  sum  of  the 
coordinates  of  a  point  was  even,  an  atom  was  placed  there,  otherwise  a  space 
was  left.   Truncation  of  quotients  in  fixed  point  arithmetic  allows  a  simple 
test  for  odd  or  even  numbers.   The  fixed  point  coordinates  of  the  atoms  were 
stored  in  memory  as  floating  point  numbers  for  later  calculations.  The 
"core"  atoms  were  assigned  numbers  from  2-64  and  the  fixed  atoms  were  number- 
ed from  65  to  172. 

In  the  calculation  of  the  forces  or  potentials  associated  with  an  atom, 
the  atoms  closer  than  nearest  neighbors  in  an  undisturbed  lattice  must  be 
considered.   A  "nested"  Do  loop  was  used  for  these  calculations.   The  out- 
side Do  loop  was  indexed  from  1  to  64  and  the  inner  Do  loops  initial  index 
was  always  one  greater  than  the  outer  Do  and  ended  at  172.  This  arrange- 
ment avoids  the  calculation  of  forces  or  potentials  twice  for  the  same  pair 
of  atoms.   The  x,  y  and  z  coordinate  separations,  in  that  order,  between 
the  atom  under  consideration  (index  of  outer  loop)  and  any  other  atom  were 
tested.   If  any  of  these  separations  were  found  to  be  greater  than  ROE, 

45 


the  inner  Do  index  was  advanced  by  one.   If  the  coordinates  passed  these 
three  tests,  the  square  of  the  distance  between  atoms  was  compared  to  the 
square  of  ROE  (R0E2).   The  true  distance  of  any  atom  that  passed  this  test 
was  used  in  the  force  or  potential  function  and  the  result  applied  to  both 
atoms.   (In  potential  calculations,  half  the  potential  was  given  to  each 
atom).   If  the  square  of  the  distance  between  atoms  was  equal  to  or  greater 
than  R0E2,  the  inner  Do  index  was  advanced  by  one.   Most  of  the  atoms  will 
be  eliminated  in  the  first  three  tests.   The  procedure  avoids  time  consum- 
ing square  root  calculations  by  the  computer. 

The  double  iteration  method  used  to  move  atoms  requires  two  calcula- 
tions of  resultant  forces  for  each  "time  step".   The  same  section  of  the  program 
was  used  for  both  force  calculations,  but  the  equations  required  to  move  the 
atoms  twice  in  the  time  step  are  not  identical.   The  fixed  point  variable 
INDEX  was  used  as  a  switching  device  to  jump  over  the  motion  equations  that 
locate  the  atoms  at  temporary  positions.  The  jump  is  executed  after  the 
second  calculation  of  forces. 

The  movement  of  atoms  is  continued  until  the  kinetic  energy  of  the 
target  is  a  maximum.   It  is  necessary  to  advance  one  time  step  beyond  the 
maximum  to  establish  the  maximum.   The  particle  locations  and  velocities  at 
the  time  of  the  maximum  are  still  in  memory,  and  these  are  used  in  kinetic 
and  potential  energy  calculations  or  as  output. 

Although  the  two  body  problem  is  in  two  dimensions,  the  impact  para- 
meter is  the  same  as  that  in  the  lattice.  The  two  body  problem  is  solved 
in  the  same  manner  as  the  lattice  except  that  the  interaction  is  stopped 
when  the  atoms  are  separated  by  a  distance  equal  to  ROE  which  means  the 
potential  and  force  between  atoms  is  zero.   The  interaction  is  initially 
started  with  the  separation  of  atoms  slightly  less  than  ROE.   The  AT  for 
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the  two  body  problem  is  much  smaller  than  the  one  used  for  lattice  calcula- 
tions.  The  small  Z^  T  consistantly  produces  very  much  more  accurate  re- 
sults in  the  two  body  problem  than  in  the  lattice. 

The  input  consists  of  the  bullet  energy  and  mass,  its  x  coordinate, 
and  the  time  step  multiplier.   In  order  to  calculate  data  for  all  36  impact 
points,  six  input  cards  are  required;  each  is  identical  except  for  the  x 
coordinate  of  the  bullet.   The  x  coordinates  assigned  are  3.0,  2.9,  2.88 
2.7,  2.6,  and  2.5  in  that  order.   After  receiving  the  input  from  one  card9 
the  corresponding  z  coordinate  on  the  line  AB  is  calculated  (see  Fig.  6). 
The  computer  proceeds  through  the  program  for  this  point,  then  subtracts  0.1 
lattice  units  from  the  z  coordinate  and  repeats.   It  procedes  down  a  column 
until  all  calculations  for  impact  points  on  the  line  BC  have  been  made.  The 
computer  pauses  at  this  point  to  await  instructions.   Here  it  may  receive  a 
new  data  card  and  start  another  column.  After  the  calculations  for  the  point 
on  the  apex  of  the  triangle  (x*2.5,  z«3.5)  are  complete,  the  computer  can 
be  directed  to  print  the  output  in  the  triangular  form  shown  by  Fig.  6. 
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APPENDIX  n 


BLOCK  DIAGRAM 


Start 


Dimension  and  common 

statements 
Eroded  force  and  potential 

subtraction  terms 
BZ  and  stopping  point 
Original  location  of  bullet 
Radius  of  effect  (ROE) 


BLOCK      A 


HD 


6 


*     \lnput ) 


Velocity  of  bullet 

Time  step 

Constants  used  in  velocity 

and  movement  equations 
Sat  lattice  velocities  *0 


BLOCK   B 


Location  of  lattice 
BLOCK   C 


Set  potential  and  kinetic  energy 

of  lattice  atoms =0 

Store  original  nositions  of  all 

atoms 

BLOCK   D 


To  BLOCK  E 
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Prom  BLOCK   D 


Set  all  forces  -  0 
Calculate  resultant  force 
on  all  atoms 


BLOCK   S 


Movement  of  atoms 
to  temporary 
positions 


BLOCK  PI 


Movement  of  atoms  to 
new  positions 

BLOCK  P2 


€t 


Calculation  of  individual 
and  total  kinetic  and 
potential  energy 

BLOCK  H 


Output 
s««  Fi&°    3  block  j 


To  BLOCK   K 
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From  BLOCK  J 


Jl. 


Location  of  atoms  for  two 

body  problem 
Calculate  velocity  of 

bullet 
Set  force  =  0 

BLOCK   K 


& 


Calculate  forces 
between  atoms 

BLOCK  L2 


Move  atoms  to 
temporary 
locations 

BLOCK  Ml 


Move  atoms  to 
new  locations 


BLOCK  M2 


C 


Output 
See  Fig. 


> 


Calculation  of  ratios  and 
parameters  to  comoare 
two  body  oroblem  with 
lattice  problem 


BLOCK   N 


To  BLOCK   R 
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From  BLOCK  N 


Storage  of  comparison 
data  for  each  impact 
point 


BLOCK  R 


<T) 


Output 
Impact  area  in 
triangular  and  tabular 

form     BLOCK  W 


End 


A 


Potential 
function 

BLOCK  Y 


E 
L2 


Force 
function 

BLOCK  Z 
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APPENDIX  IE 


DETAILED  FLOW  CHART 


START 


I 


DIMENSION 
AND 
COMMON 

STflfTEMfchTTS 


BLOCK  A 


save  memory  locations  for  subscripted 
variables 


L«\ 


ROE2*2.0 


impact  point  number 

lattice  unit  distance  between  nearest  neighbors 
and  the  square  of  that  distance 


EPST-  FORqftfc) 


FORCE 


ePSTP^fbTF(ROe) 


^6teiOTiAL 


3         IN^T 
E\/       BX,  DTI,  PMASS 


5    PXO)=BX 
PfcE&fc- 100,0 

cycLE  =  o.o 

INDEX-  O 

ryO)-  -^5 


I 


TO 

BLOCK  B 


terms  subtracted  from 
calculated  potential  and 
force  that  "erode"  them 


starts  here  for  every  "c 
in  impact  area 


starts  here  for  every  impact 
point;  assigns  coordinates 

to  bullet 
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V*0)=  o.o 


v,Y(,|.\  EV<faU54    *k*\\*\A 


[000.OH  PMftss 
\te(0  =  O.O 


fyi\a  \0S.  -AloBqU   v   lO"27 
DTOC=    DT/3.61M  x  \0"'° 
DTORN\=    CJT/R^ 
0TO2RM  -   OTORN\/2,0 
RMt=.PW^*  l.bS^^JO 
PTORMJ  r   DT/RMI 
PT02RAM  aJ*TORAU/fe.O 


-10 


r27 


VX$)=  O.O 


£50 
COMT\KJUE 


TO 
BLOCK  C 


BLOCK  B 

velocity  components  of 
bullet 


calculates  constants  used 
in  force  and  movement 
equations 


sets  velocity  of  lattice  atoms 
to  zero 
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(si  ^b5 


E 


T-- 


*_ 


DO  hO  J-  Ijl 


BLOCK  C 
first  numbers  to  be  assigned  to  lattice  atoms 
2  for  movable  core,  65  for  fixed  atoms 


uy 


IT,  JT,  KT  are  X,Y,  Z  coordinates 
in  the  system  used  for  reference 


^cr=o 


DQ5q  k*  \njf^@ 


lT=o 


^_ 


DOS9  Is  1/7 


<— (S8 


NTC=  IT*  JT  t  KT 


each  point  with  integer  values  for 
coordinates  is  tested;  if  the  sum  of 
the  coordinates  is  an  even  number, 
an  atom  is  assigned  to  that  position 


these  decision  processes 
determine  if  an  atom  is  in 
"core";  this  is  necessary  in 
order  to  have  core  atoms  with 
consecutive  numbers  2-64 
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±o 


JT 


<o 


Zo 


:jt-^ 


>0 


60 


M=M  + 


end  of  Do  loop  over 
Y  coordinates 


30 

RK(N)  =  IT 


tf«tf-H 


£. 


COKTTlMUE 


KT*tCT  +  \ 


5^ 

CONTINUE 


I 


JT=JT+ 


I0O 

COKJTl»JUt 


BLOCK  C  (  con't) 


for  a  given  Y  and  Z,  X 
coordinates  are  checked 
first;  then  Z  is  changed 
X  coordinates  checked  again; 
after  tests  on  X  and  Z  , 
Y  coordinate  is  changed 


end  of  Do  loop  over  X  coordinates 


end  of  Do  loop  over  Z  coordinates 


t6 

BLOCK  "D 
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/C 


fPE(l)-  O.O 
R1K(I)=R*(D 

RlY(l)=RY(I) 
WB(I)--R2(I) 


US" 
CowtikJuE 


^[pO  UU    I  *    I,  (oH 

?PE(I)-  OiO 
PKE  (I)  =  O.O 

RXR(I)-  R*(l) 
RYR(I)=RY(J) 


R?R(i)=R2a) 


TO 
BLOCK  E 


BLOCK  D 


sets  temporary  and  permanent 
coordinates  of  immovable  atoms 
equal  to  eachother 


sets  initial  potential  and  kinetic  energy 
of  core  atoms  equal  to  zero 


stores  original  positions  of  bullet  and 
core  atoms 
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70 

F?(l)  =  O.O 
I 


BLOCK  E 


initial  forces  between  core  atoms 
set  to  zero 


15" 

COKTH^OE 


do  aoo  i--  l,W 


IP1=>\ 


Do20oJ«IPl/l2 


4^I-\ 


<o 


^0 


>o 


keeps  bullet  from  interacting 
with  immovable  atoms  covering 
bombarded  face 


D*X-R*(J)-SX(J) 


ftt)RKsftBSF(DfcK) 


tests  X  coordinate 
separation  between  atom 
I  and  J 


io 


20O 


<o 
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<o 


PRY=RV(J)-RY(I) 


rdry  =  ^bsf(dky) 


"^DRY-RO^V- ^(J 


i     ; 


D1S7=DR*2+  D*YV  DR21* 


DIST-AJPIST 


BLOCK  E  (con't) 


test  Y  coordinate  separation 
between  atom  I  and  J 


test  Z  coordinate  separation 
between  atom  I  and  J 


test  square  of  distance  between 
atoms;  if  square  passes  test, 
actual  distaice  is  found  and  used 
to  calculate  force 


FORCk*T=ORC(QIST)-£PST  « 


FORCE 
"FU  OCT  |  OK) 


FOD:  FOfcC£/D3ST 
FR  =  PRX*  FOD 
FB  =  "DRY  *  FOD 
FC=DR2   *  FOD 


reduce  total  force  into  X,  Y,  and  Z 
components 
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BLOCK  E  (con't) 


\J-fe5 


-CO 


>o 


fx(j)=f*Ct)i-fa 
4fy(j)=py(j)+f& 


Fk(i)  =  fx(j;-fa 

FYU)a^(l)-FB    \t 


add  force  to  total  forces  on  atoms 
I  and  J  if  they  are  in  the  core 


-zoo 


end  of  Do  loop  for  finding  forces 


r->D02>3S*I  =  l>M 


>0  * 


>0  * 


2ZS 

COVJTIK)UE 


So 


ypcttjs  0.0 


*Pf(l)*O,0 


^S<f^(J))>l^F2(i)^0 


TO 

Block  Fl 


round  off  procedures  in 
computer  may  lead  to  very 
small  resultant  forces  in 
situations  where  the  forces 
should  be  zero;  this  procedure 
sets  to  zero  any  resultant  force 
less  than  10       newtons 
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BLOCK  Fl 


DO  300    1-2,  <©4 

Rixd)  =  i^x^cx; 

FH(I)--  FA  (I) 
Vlx(l)  =  M*tt) 


store  position,  velocity,  and  force 
of  core  atoms 


VK(I)-FX(I)  *WORM  *  M*(I) 
a\s>o  SQW.C  cslc.  -foir   v^^  ?  eonrtporKMk 


calculate  temporary  velocities 


move  core  atoms  to  temporary 
locations 


3oO 
COvmMUE 


«Uo  same,  pirocedo^e 
•for   S  £  '  COiwpo^cmrJs 


1 


store  position,  velocity,  and  force 
of  bullet 
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BLOCK  Fl  (con't) 


VX(i)=  FXCi)^DTORMli-MK(l) 
Also  same  -Pov  L|  £  ^.  corwpo^eMs 


also   scurvy  -for   u  £^    coo^*/>q4cs 


in1dex=1 


SOTO   10 

(in  block e) 


since  mass  of  bullet  can  be 
varied,  equations  of  motion 
for  the  bullet  must  be 
calculated  separately 


recycles  to  calculate  forces  at 
temporary  locations 
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FROM  IfODBX  check 
Block    F-l 


r* 


BLOCK  F2 

forces  at  temporary  locations  have 
been  calculated 


DO  +25     I  -  2y  1*4 

RX(T)=[\)K(I)  +  VJU(1)]^DT0C  iRl\(l) 
ak>o  SQry^    calculation*  -Co*-  u  £  ^   cowponetvls 


4?s 

Cov)Tioue 


I 


IWDFX^O 


CYCLE*  CYCL£+\.0 


TO 
3LOCK.  S 


new  velocities  calculated 
using  average  forces 

new  positions  assigned 
using  average  velocities 


separate  bullet 
calculations 
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PKE8*32.SlUM06»[vx(ef  +  MYtef+M^f] 


BLOCK  G 
kinetic  energy  of  target 


1 

ho 

PKE8fc=PKE8 

<L                             y* 

> 

t 

GO  TO    "70 

(JN  Block  e) 

<o 


if  kinetic  energy  of  target  greater  than 
former  value,  recycle  to  calculate  forces 
and  move  atoms  again 


kinetic  energy  of 
target  compared  to 
former  value 


if  kinetic  energy  is 
less  than  former 
value,  movement 
of  atoms  is  stopped 


Covrn>JUE 


TO 

B^ocx  H 
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DI*  =  RiX(9)-MX&) 
DlY=ttY(6)-RlY&) 


BLOCK  H 


coordinate  separation  between  target  and  bullet 


D81=D1X*hDIY2+DIZ* 


distance  between  bullet  and  target 

is  found  and  potential  energy    calculated 


^D8l^^81 


poTl&=o.o 


±L 


P0Tl^cpcrrF(D&l)-EPSTP^ 


PoTcun^L 
po»JcT«on 


PKEl^  0.5i8  *  \d8*[\|l*(l)%  viY(l)\  VJZfl)*]*  PMftSS 


kinetic  energy  of  bullet 


Pte(l)--PKE1 
T?OT=  O.O 
JPK£=  O.O 


DO  bOO  1^  !>4 


1P1-I+1 


calculate  potential  energy  of  each  atom 
and  the  total  potential  energy 


DOboo  j*ipi,n2 


i-l>^%s^o 

<c 
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BLOCK  H  (con't) 


SO 


>-^((oOO 


prevents  interaction  of  buLet  with 
immovable  atoms  covering  bombarded 
face 


^L^K=RU(J)-R1X(I) 


adrx-absf(dr>) 


DRY=R1Y(J)-R1Y(I) 


APRY=A&SF(DZY) 


DR2  =  R1?(J)-  R1Z(I) 


fiDRt'-BGSFiVRl) 


to 


boo 


<o 


test  X  coordinate  separation 
between  atom  I  and  J 


test  Y  coordinate  separation 
between  atom  I  and  J 


test  Z  coordinate  separation 
between  atom  I  and  J 
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<o 


DJST=  DRX2+  DRY  V  D'Rfc'1' 


^O 


BLOCK  H  (con't) 

square  of  distance  between  atom 
I  and  J 

test  square  of  distance,  if  it  passes 
test,  find  actual  distance  and  ase  it 
to  calculate  potential  between  I  and  J 


<o 


DIST^DIST 


I 


POT=POTF(DIST)  -  EPST? 


PPE(I)=  PPE(I)-KD.5*POT 
P?E(J)=  PP£  (J)  +0.5*  POT 


POTEMTlAL 

Fuuctiokj 


half  of  potential  between  atoms  is 
assigned  to  each  atom  (  I  and  J  ) 


TFOT=  TPoT  +  PtfT 


potential  energy  summation 


WOO 
cotomoue 


end  of  Do  loop  for  calculating  potential  energy 


COL3=  pPEdJ-O.S'^POTlS 
pE8A  =  PPE(&)-0.5*Pcnrl8 


I 


DO   bSO      1=  2,  feV  t  . 

PKE(l)»32.^lbMO-8[VlX(l)  4  VJY(I)  +Vl*(l)  J 

TPKE  =  TPKE  +  PKE(I)  


kinetic  energy  of  core 
atoms  not  including  the 
bullet 


Coajt/aJoe 


1 

TPKE=TP<E  +  PKE1 
TF^TPKE-hTPOT 

\t 

I         EV 

total  kinetic  energy  and  total  energy 


percent  error  in  total  energy 


TO 

Block  J 
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1 


BLOCK  J 


OUTPUT 
see    p;3.   4 


BLOCK K 


COKiV  =•    O.Ol/DTX 
PTOC  =  PTOC  *  COUV 

PTORNA  s  PTORM  *COkJV 

PT02QNA  -  DT021M*  CONV 

ptoRmi-dtoRmuccomv 

PT02  fc  Ml  --PTQ2RM3  *COkJV 


change  constants  for  two  body  problem; 
DTI  always  equal  to  0.01  in  two  body 
problem 


RK2-  O.O 
RY2--0.0 


target  located  at  origin 


DSQ  =  (b*  -3.0  f  +(BX-  3. 0) 


square  of  impact  parameter  in  the  lattice 


PXl  =  ^RO£2-DSd    -  10** 
I 


VX2-*  O.O 
VV2s  O.O 

vyi=  o.o 
pxl*  o.o 


FX2-  O.O 

FY1--  o.o 
FY2  ■  o.o 


assign  bullet  Y  coordinate  to  give  proper 
impact  parameter 


calculate  X  coordinate  of  target  to  insure 
separation  of  bullet  and  target  is  slightly  less 
than  nearest  neighbors  separation 


itooexrO 


TO 

Block  l1 


assign  velocity  to  bullet 
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BLOCK  LI 


\OOS 

DX  =  RX1-  RX2 
DY  =  RY1-RY2 


DI5T=  DX'  +  DY' 


calculate  distance  between  atoms 


>o 


\CHO\       test  to  see  if  less  than  ROE 


PIST=  NlOIST 


Force-  fo  rc(-dist)-ePst^ 


Riioctiok) 


BLOCK  L2 


calculate  force  between 
atoms 


P*l=  (dx/pist)*  FORCE 

FX2*  -FM 

FY1 a  (pv/D\%T)*Fofece 
FY2=-FY1 


total  force  is  divided  into  X  and  Y 
components 


Block  m£ 


68 


RX1T-    RX1 
RX2T  =    RX2. 

FXlT  £  FXl 
VX1T  =  VX1 
VX2T  =  VAZ 
cx\%o  s>a*v\e  -Po*-   u 


COi 


BLOCK  Ml 


store  old  positions,  forces,  and 
velocities 


VX1  -   FXL*DTORMl    -hVXl 
V*2r  1pX2  *  DTOPm     4-  VK2. 
Rxl  -  QvJxl  +  N/XlT]  H  PTOC  +-  rx,it 
K*2-  [VX2  VVX2TJ  *  DTOC  +  RV2T 

Q\so      Y    Go^po^ca4s 


ind^x*  i 


GO  TO 
10O5" 

Block  Ll 


FROM  i wpbk  test -Block  L2 


i 


VXl   *  ^FXl  +  FXlT)*DT02ftMl  +  \JXlT 
N/XZ  =  (FX2  +  FK2T)*pT02RKi-\M2T 
RX1  =  (\JX1  4  VX1T)  *  DTOC  +■  ftXlT 
K*2.  =  (vJX2  +  VX2T)  *  DTO  C  +  £x  2T 

olso    V    Co/*npovoef>"Ks 

\ 

r 

1 1 

lUDcy  —  ^v 

^» 

&o   TO 
IOOS 

&LOCX  Ll 

\Nver 

-\  -  w 

^ 

calculate  temporary 
velocities  and  positions 


recycles  to  calculate  forces 
at  temporary  positions 


BLOCK  M2 


calculate  final  velocities 
and  positions 

average  forces  and  velocities 
are  used 


recycles  to  distance  test 
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J07D 

Test    SLoc/r/_l 


I 


BLOCK  N 


TKE  -  ?<?KE.l4-P2KEa 


Out  Pu 


CALCuu^-TlOM    OF  SeOeAAl 

VALUES  AfOO  fcATlOS 
USE"D  TO    COM?P»P^ 
LftTTICg     U){TH  TWO  feODY 


c 


OUT  PUT 
See    "PiG.    5 


kinetic  energy  of  target  and 
bullet 

total  energy  of  system 


TO 

Block  R 


BLOCK R 


STORE    \/ftftlOU&    VALOe* 
fOO/JP    FOC    LATTICE:    fl»JD 
TV^>0    fcO£>y    "PROBLEM     i*JTu>0 

Pi  toe /usionJal.    ApR^ytcftO] 


index  I  is  value  stored 

index  L  is  impact  point  number 


impact  point  number  increased  by  one 
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B2^B2-o.l 


COKJTI  K)Ug 


BLOCK  T 

Z  coordina  te  of  impact  point  lowered  by 
0   lOlattice  units 


Z  coordinate  of  impact  point 
tested  to  see  if  it  is  the  end 
of  a  column 


TO 

Block  Uj 


if  sense  switch  one  is  up,  new 
data  card  is  read  to  start  another 
column  in  "triangle" 


BLOCK    W 


1 


O  UT  PUT 

OF   OftTft 
see   Pis.  (o 
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Function  Potf(d\st) 

CON\rAOKJ 

£TFP,  APfiC  fTFF 


BLOCK  Y 


-go^H  p3-^»H 


Potp  ■  (etff)(i.Q4-[pIst^prcT  6°lH) 


I 


RetorO 


JL 


e  k)d 


,1* 


PUKJCTIOK)     FORC(DIST) 
COIAMOkJ 
eTF-P,  AFAC;  tTFF 

>' 

,.-'                                                             L 

yp=l.o+Y 

I 


?s(yp) 


-3.131 


BLOCK  Z 


Pottc»fcTrp)«(z)«(i.o  +  arw^J 


OlST-*  PIST 


I 


RetoUU 

I 


1  EMD 
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APPENDIX  IV 
DEFINITION  OF  VARIABLES 

RX(I)  X  coordinate  of  atom  number  I  at  any  time  (lattice  units) 

RlX(I)         x  coordinate  of  atom  number  I  at  former  position  (same  as 
position  1  in  procedure  section),  stored  wbile  average 
forces  are  computed  (lattice  units) 

RY(I),  RZ(I),  R1Y(I),  R1Z(I)   y  and  z  coordinates  defined  in  the  same 
manner  as  RX(I)  and  R1X(I) 

FX(I)  x  component  of  force  on  atom  number  I  at  any  time  (newtons) 

FlX(I)         x  component  of  force  at  former  position  of  atom  number  I 
(newtons) 

FY(I),  FZ(I),  F1Y(I),  F1Z(I)   y  and  z  components  of  force  defined  in  the 
same  manner  as  FX(I)  and  FlX(I) 

VX(I)  x  component  of  velocity  at  any  time  of  atom  number  I  (m/sec) 

V1X(I)         x  component  of  velocity  of  atom  number  I  at  former  position 
(m/sec) 

VY(I),  VZ(I),  V1Y(I),  V1Z(I)   y  and  z  velocity  components  defined  in  the 
same  manner  as  VX(I)  and  V1X(I) 

C(I,L)         matrix  for  storing  data  I  from  impact  point  L 

PPE(I)         potential  energy  of  atom  number  I  (ev) 

PKE(I)        kinetic  energy  of  atom  number  I  (ev) 

RXR(I),  RYR(I) ,  RZR(I)    x,  y,  and  z  coordinates  of  all  atoms  in  initial 
positions  (lattice  units) 

ETFP,  AFAC,  ETFF  constants  used  in  Thomas-Fermi-Firsov  potential 

L  impact  point  number 

ROE  nearest  neighbors  separation  (lattice  units) 

ROE 2  ROE  squared 

EPST  value  subtracted  from  calculated  force  to  erode  the  force 

(newtons) 

EPSTP         value  subtracted  from  calculated  potential  to  erode  the 
potential  (ev) 
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FORC 

POTF 

EV 

BX 

DTI 

PMASS 

BZ 

BZS 

PKE8R 

CYCLE 

INDEX 

DT 

RM 

RM1 

DTOC 


DTORM 

DT02RM 
DTORM1 

DT02RM1 

M,N 

JT 


function  subroutine  that  calculates  force;  argument  is 
the  separation  of  atoms  (newtons) 

function  subroutine  that  calculates  potential;  argument  is 
the  separation  of  atoms  (ev) 

original  kinetic  energy  of  bullet  (ev) 

original  x  coordinate  of  bullet  (lattice  units) 

time  step  multiplier 

mass  of  the  bullet  (amu) 

original  z  coordinate  of  bullet  (lattice  units) 

z  coordinate  of  bullet  at  end  of  column  in  impact  triangle, 
used  as  stopping  point  (lattice  units) 

kinetic  energy  of  target  at  the  end  of  previous  time  step 
(ev) 

number  of  time  steps 

switching  device 

length  of  basic  time  step  (seconds) 

mass  of  copper  atom  (kg) 

mass  of  bullet  atom  (kg) 

a  constant  used  in  movement  equations;  when  multiplied  by 
the  sum  of  two  velocities,  it  computes  the  average,  multiplies 
it  by  DT  and  gives  the  result  in  lattice  units  (sec  x  lattice 
units/m) 

a  constant  used  in  velocity  equations;  it  combines  the 
quotient  DT/RM  into  one  constant  to  avoid  repeated  divisions 
(sec/kg) 

same  as  DTORM  except  it  averages  the  force 

same  as  DTORM  but  used  in  bullet  calculations  only,  since 
the  bullet  mass  may  be  different  from  that  of  the  lattice 
atoms 

same  as  DT02RM  but  used  in  bullet  calculations  only  since 
bullet  mass  may  be  different  from  that  of  the  lattice  atoms 

atom  numbers 

y  coordinate  of  an  x-z  plane 
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KT  z  coordinate  of  an  x-y  plane 

IT  x  coordinate  of  a  y-z  plane 

NTC  sum  of  x,y,  and  z  coordinates  of  a  point  in  the  lattice 

DRX,DRY,DRZ    difference  between  x^y,  and  z  coordinates  respectively  of 
any  two  atoms  (lattice  units) 

ADRX,  ADRY,  ADRZ  absolute  value  of  DRX,  DRY,  and  DRZ 

DIST  distance  between  two  atoms  or  the  square  of  the  distance 

(lattice  units) 

eroded  force  between  two  atoms  (newtons) 

FORCE  divided  by  the  distance  between  two  atoms  (newtons/ 
lattice  unit) 

x,y,  and  z  components  respectively  of  the  force  between  two 
atoms  (newtons) 

kinetic  energy  of  target  atom  in  the  lattice  (ev) 

x,y,  and  z  coordinate  separation  between  target  and  bullet 
(lattice  units) 

distance  between  target  and  bullet  in  lattice  or  the  square 
of  the  distance  (lattice  units) 

potential  between  target  and  bullet  (ev) 

kinetic  energy  of  bullet  at  end  of  time  step  (ev) 

total  potential  energy  of  lattice  (ev) 

total  kinetic  energy  of  lattice  (ev) 

potential  between  any  two  atoms  (ev) 

potential  of  the  bullet  with  respect  to  the  lattice  minus 
the  target  (ev) 

potential  energy  of  the  target  with  respect  to  the  lattice 
minus  the  bullet  (ev) 

total  energy  of  the  lattice  (ev) 

percent  error  in  total  energy  (percent) 

time  step  multiplier  conversion  factor  for  two  body  problem; 
makes  two  body  time  step  multiplier  equal  to  0*01  regardless 
of  the  value  used  for  the  lattice  problem 


FORCE 
FOD 

FA,FB,FC 

PKE8 
DIX,DIY,DIZ 

D81 

POT  18 

PKE1 

TPOT 

TPKE 

POT 

C0L3 

PE8A 

TE 

PERE 

CONV 
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DSQ  square  of  the  impact  parameter  in  the  two  body  problem 

(lattice  units  squared) 

XIP  impact  parameter  in  the  two  body  problem  (lattice  units) 

RX2,RY2        x  and  y  coordinates  of  target  at  any  time  in  two  body 
problem  (lattice  units) 

RX2T,  RY2T      x  and  y  coordinates  of  target  at  former  positions,  stored 
while  average  force  is  computed  (lattice  units) 

RX1,RY1,RX1T,RY1T   x  and  y  coordinates  of  the  bullet  in  the  two  body 
problem  defined  in  the  same  manner  as  RX2SRY2S RX2TS  and 
RY2T 

VX2,VY2        x  and  y  components  of  target  velocity  in  two  body  problem 
(m/sec) 

VX2T,vT2T      x  and  y  components  of  target  velocity  in  two  body  problem 
at  former  position  (m/sec) 

VX1,VY1,VX1T,VY1T  x  and  y  components  of  velocity  of  the  bullet  in  the  two 
body  problem  defined  in  the  same  manner  as  VX2,VY2SVX2T9 
and  VY2T 

FX1,FY1        x  and  y  components  of  force  on  the  bullet  in  the  two  body 
problem  (newtons) 

FX1T,FY1T      x  and  y  components  of  force  on  the  bullet  in  the  two  body 
problem  at  its  former  position  (newtons) 

FX2,FY2,FX2T,FY2T  x  and  y  components  of  force  on  the  target  in  the  two 
body  problem  defined  in  the  same  manner  as  FX18  FYljFXlT, 
and  FY1T 

DX  x  coordinate  separation  between  atoms  in  the  two  body  problem 

(lattice  units) 

DY  y  coordinate  separation  between  atoms  in  the  two  body  problem 

(lattice  units) 

P2KE1         kinetic  energy  of  the  bullet  in  the  two  body  problem  after 
interaction  is  complete  (ev) 

P2KE2         kinetic  energy  of  the  target  in  the  two  body  problem  after 
interaction  is  complete  (ev) 

TKE  total  kinetic  energy  in  the  two  body  problem  (ev) 

T2RA  recoil  angle  in  the  two  body  problem  (degrees) 

B2SA  scattering  angle  in  the  two  body  problem  (degrees) 
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SRE8 

SRE1 

ATLRA 

TLRA 

ABLSA 

BLSA 

RBMT 

EABL 


resultant  velocity  of  target  in  lattice  (m/sec) 

resultant  velocity  of  bullet  in  lattice  (m/sec) 

quotient  of  two  velocities;  used  to  find  TLRA 

recoil  angle  in  the  lattice  (degrees) 

quotient  of  two  velocities;  used  to  find  BLSA 

scattering  angle  in  the  lattice  (degrees) 

bullet  to  target  mass  ratio 

energy  absorbed  by  the  lattice;  does  not  include  potential 
or  kinetic  energy  of  target  and  bullet  (ev) 


ET2,ETL,RET,RKEP,RSA2L,RRA2L,RBFI  ratios  used  to  compare  lattice  interaction 
with  two  body  problem;  the  ratios  are  in  terms  of  variables 
already  defined. 
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APPENDIX 


PROGRAM  LISTING 


BLOCK  A 

PROGRAM  TFF 

CIMENSION  RX(200)T  RY(20C),  RZ(2CC) 

DIMENSION   FX(2CO),FY(2CC).:  Z(20C) 

DIMtKSIUN    VX12CO)  ,  VY(2CC)  ,  VZ(2C0) 
OOIMENSION   R1X(2C0),  RlY{2CC),  R1ZC200),  FlX(20C),  F1Y(2CC), 
1   F1Z(200)  ,  VlX(200),  V1Y1200),  V1Z(2CC) 

DIMENSION  C( 22,36  )  ,PPE(2CO) ,  PKE ( 200 ) , RXR (200 ) , RYR ( 200 ) , PZR ( 200 ) 

COMMON  ETFP,  AFACETFF 

AFAC=3.5906 

ETFF=5.9360E-6 

ETFP=6.6980E+3 

L=l 

ROE  =  SQRTF  (2.C) 

R0E2  =  2.C 

EPST  =  FORC  (  ROE  ) 

EPSTP    =    POTF(ROE) 
10       F0RMAT(34X5UH    THOMAS-FERM I-F IRSCV    POTENTIAL,     ERCCED,    COPPER-COPPER 
1  ) 

3    READ    INPUT    TAPE    2,4,    EV,8X,DTI,    PMAS S 
4       F0RMATUE10.3) 

BZ  =  BX  ♦  1.0 

BZS  =  6.0-BX 
5    RX( 1 )  =  BX 

RZ( 1)  =  PZ 

PKE6R  =  -100.C 

CYCLE  =  0.0 

INDEX  =  0 

RY( 1 )  =  -SQRTF(2.0) 
12   0F0RMAT(/,29H  BULLET  INCOMING  LOCATION   X=  F6.3,  5H    Z=  F6.3, 

110H    ENERGY=  F8.1,  2UH    TIME  STEP  MULTIPLIER=   F5.3,  /   ) 
15    FORMAT(lHl) 

BLOCK  B 

V  X  (  1)  =  0  .  c 

VY<1)=  SQRTF(  EV  »  63 . 54 / ( 1 . CE2  *  PMASS)  )  »  5.511811  EOH 
i-  VZ(  1  )  =  0.0 

DT  =   DTI/  VYM)  »  1.807E-1C 
RM=  105. 463911  E-27 
DTOC   =  DT/3.6H4  E-10 
DTORM  =   DT  /  RM 
DT02RM  =  DTORM/2.0 
RM1  -    PMASS  »1.6598  E-27 
DT0RM1  =  DT/RM1 
DT02RM1   =  DT0RM1/  2.0 
DO  250  I  =  2,  64 
VX{  I)  =  0.0 
VY(  I  )  =  0.0 
250  CONTINUE 
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PLCCK  C 


20 

2 

2 

2 

2 

25 

50  RX(N) 

RY(N) 

RZ(N) 

N  =  N 

GC  TO 

rx ( w )  =  i  r 

R  Y  (  M  )  =  J  T 
RZ(M)  =  KT 

m  =  M  +  l 

I T= I T  +  1 

CONTINUE 

KT  =  KM 
CONTINUE 
JT  =  JT  + 
CONTINUE 


40 


57 
58 


59 


60 


+  I 
57 


1 


RLOCK  0 


DO  6  5  I  =  65i  172 
PPE(  I  )  =  0.C 
Rl X(  I  )  =  RX(  I  ) 
R1Y(  I )  =  RY(  I  ) 
R12(  I  )  =  RZU  ) 
6  5   CONTIMJt 

CO  66   1=  1  ,6*1 
PPE<  I  )  =  0-C 
PKE(  I  )  =  C.C 
RXRK  I  )  =RX( I  ) 
RYR(  I  )  =  RY (  I  ) 
66  R2R( I  )  =  RZ( I ) 
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PLCCK  E 


70 


75 


130 
135 

120 


mo 

1U5 

150 
155 


160 
165 

200 


21C 
215 
220 
225 
230 
235 


IP 

on 

IF 

IF 


DC  7S 
FX(  I  ) 
FY(  I  ) 
FZ  (  I) 

cnNTi 

00  20 

1  = 

20 

( I 
<  J- 

1F{  J- 
DRX  = 
ADRX 
IF  (AD 
DRY  = 
ADRY 
IF(  AL; 
DRZ  = 
ADRZ 
IF(  AO 
CIST 
IF  ( 

01  ST 
FORCE 
FOO  = 
FA  = 
FR  = 
FC  = 
LF( 
FX( 
FY( 
FZ(  J) 
FX{  I  ) 
FY( 
FZ( 

cw 

DO 

IF( 
FX( 
IF( 
FY{ 

IF( 
FZ(  I) 
C0N1TI 
IF  (  I 


I 

=  C 


=  1 

« c 

o.c 

0.0 


6U 


J 
J) 
J) 


I  ) 
I  ) 
TI 
23 

A 
I) 

A 
I  ) 


MJL 
0  I 

I 
C  J 
-1  ) 
65) 
58) 

PX 
=  A 
PX 

RY 
=  AH 
RY  - 

RZ( 
=  Afl 
RZ  - 
=  (I 
DISI 
:  SL 

=  F 

FOK 
DRX 
DRY 
DRZ  * 
-6  5)  1 

=  i;  X( 


=  1  ,fb 
1 

U 

}? 

?  ^ 

J) 

SF 

R 

J) 

SF 

r> 
J)' 

SF 

R 
RX 

RT 
CR 

CE 

6 


=  F  Y( 


NUE 
5  I 
H^F 

H.SF 

ABS 

NUE 

NCE 


IP 

0, 
0, 

r\ 
^  T 

( 

OE 

(I) 

cr 

(0 

OF 
* 
p 

F( 

C( 

/ 

FC 
FC 
FC 
60 
J) 
J) 
J) 
I  ) 
I  ) 
I  ) 


30, 
=55, 

200 
RX( 
RX) 

1U 
RY( 
Y  ) 

}U 
RZ( 
I) 

15 
OR* 
E2) 
CIS 
0  I  s 
OIS 


17 
12 
13 
i 

I  ) 

C 
I  ) 

5, 
I) 

C 


0 

c; 

12C 


200,2C0 

2  0  0  ,  2  C  C 


T) 
T) 

T 


200,200 

+  CRY* DRY  +   ORZ  *  DRZ  ) 

155,  200,  200 

-  r.PST 


165, 
+  F  A 
+  FC 
+  FC 
-FA 

-  FB 

-  FC 


16: 


=  1,  oU 
( F X (  I  ) )  -l.Cc-IU)  210,  210,  215 
CO 

(FY(I))  -l.Ct-lU)  220,  220,  225 
CO 

I  (FZ(  I )  )  -1.CL-1U)  230,  230,  235 
CO 

x)   2U0,  2UCt    U2C 
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FLOCK  Fl 


2UC 


30C 


DO 
Rl 

S) 

Fl 
FT 
Fl 
VI 
VI 
VI 
VX 
VY 
VZ 
RX 
RY 
RZ 

co 

Rl 
Rl 


R 

F 

F 

F 

V 

VI 

VI 

VX 

VY 

VZ 

RX 

RY 

RZ 

IN 

GO 


30  C 
X(  1  ) 

m 

X(  I 
Y(  1 
Z(  I 
X(  1 
Y(  I 
Z(  I 
(  I  )  = 
(  I)  = 
(  I  )  = 
(  I  )  = 
(I)  = 
(  I  )  = 
NT  INU 
X(  1 
Y(  1 
Z(  1 
X(l 

Yd 

Z{  1 
X(  1 
Y(  1 
Z(l 

( 1 ) 
( 1 ) 
( 1) 
{ 1) 
( 1) 
( 1) 

DEX 

Til 


I     = 


2,     6>J 
RX(  I  ) 
RY( 

fcZ 

FX 

FY 

FZ 
VX( 
VY( 
VZ  (. 
X(  I 
Y(  I 
Z(  I 

VX 

VY 

VZ 


H. 

I  ) 
I  ) 
1  ) 
I 

) 
) 


=  i- 

=  F 

=  F 

=  ( 

=  ( 

=  ( 

=  1 


RX( 
RY( 

RZ  ( 
FX{ 
FY( 
FZ( 

VX( 

VY( 

VZ( 

X{  1 

Y(  1 

Z(  1 

VX 

VY 

VZ 


*  HTCRM 

«  DTORM 
*DTCRM 

I  )  + 

I  )  * 

I  )  +    VI 

) 


1  ) 
1  ) 
1  ) 


V1X(  I)  ) 

VI Y(  I  )  ) 

(  I  )  ) 


VX(  I  ) 
VY(  I  ) 
VZ(  I  ) 

*  DTOC 

*  r:TOC 
»    DTOC 


0  TOR Ml 
DTORM 1 
DTCRM  1 
+     VlX(  1  ) 
+     V  1  Y  (  1  ) 
+     V1ZU  ) 


VX(  1  ) 
VY(  1  ) 
VZ(  1  ) 

*  DTOC 

•  DTCC 

*  f.tcc 


r  l  x  (  I  ) 
R1Y( I ) 
R1Z(  I  ) 


R  1X(  1  ) 
R1Y(  1  ) 
R1Z(  1  ) 


7C 


BLOCK    F2 


i420 


U25 


U55 

U6C 

U65 


DO     M/D 

VX(  I 
VY(  I 
VZ(  I 
RX(  I 
RY(  I 
RZ(  I 
CONT 


VX{  1 
VY(  1 
VZ(  1 
Rx(  1 
RY{  1 
RZ(  1 
INDEX 
CYCLE 


=  ( 

=  { 

=  ( 

=  ( 

=  ( 


I  = 
FX{ 
FY( 
FZ( 
VX( 
VY( 


-    (     VZ( I ) 


MJE 


(FX( 1 ) 
(  FY{  1  ) 


=  ( 

=  ( 

=  ( 

=  ( 


FZ( 
VX( 
VY( 
VZ( 


1  ) 
1  ) 
1  ) 
1  ) 


6U 

+ 

+ 
f 

+ 

+ 
+ 
+ 
f 
f 


=    CYCLE 


FIX 

n  y 

F1Z 

V  1  A 

V1Y 


(  I 
(  I 
(  I 
(I 
(  I 
(  I 


F1X( 1 ) 
F1Y(  1  ) 

F  1Z( 

VI  A 

VlY 

VIZ 

1.C 


1 
(  1 
(1 
(  1 


) 


*DTC2RM 
•DTC2RM 
»DT02RM 

*  DTOC 

*  DTCC 

*  C  TOC 

*DT02RM1 
*  DTC2RM1 

*  UT02RM1 
J  *  DTOC 

)  *CTCC 

*  DTOC 


RLOCK  G 

PKE8  =32.9163  2*8t-CP*(VX(C)*VX(8) 

IF(PKt8)  U60,  U60,U55 

IF (PKE8R-PKE8)   U6C,U65,U65 

PKE8R  =  PKE8 

GO  TO  7C 

CONTINUE 


V1X( I  ) 
V1Y(  I  ) 

V1Z(  I) 
RlX(  I  ) 
R1Y(  I  ) 
R1Z(  I  ) 


VI  X( 
V1Y( 1 ) 

V1Z(1 
R1X(  1  ) 
R1Y{  1  ) 
R1Z(  1  ) 


1  ) 


+  VY(8)*VY(8)  +  VZ(8)«VZ(8 )  ) 
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BLCCK  II 


U30 

U35 

440 


525 
53C 

5^C 


60 
I 
J 

J 


545 
550 

555 

556 

560 

600 


650 


("SIX  = 

OIY  = 

L  ■  I  Z  = 

08  1  = 

if(  nr 

D81  = 
P0T18 

GO  TO 

pot  n 

PKr  1 
VIZ  (  1 

P  K  E  (  1 

TPCT 

TPKE 

no  60 
IP  1 

CO 

If  ( 

IF( 
IF( 

CRX  = 
ACRX 
IF  (AC 
CRY  = 
AORY 
IF(  AO 
CRZ  = 
ACRZ 
I  F  (  AO 
CIST 
IF(DI 
CIST 
POT  ■ 
PPE(  I 
PPE(  J 
TPOT 
CONTI 
C0L3 
PF8A  = 
DC  65 
PKL(  I 
TPKE 
CONTI 
TPKE 
TE  = 
PERE 


R  1  A  (  8  )  -  R 1  X  (  1  ) 

R  1  Y  (  P  )     -  P  1  Y  (  1  ) 

R  1 Z  (  8  )  -  R 1  Z  (  1  ) 

(C1X*0IX  +  CIYkOIY  +  0IZ*OIZ) 
1  -  ROE" )  42  8,J*35»U35 

SQ!UFTD8l  ) 

«  PDTF  (  D8D-  -  EPSTP 

UMC 

=  CO 
=  C.5  1904  1C?L-C6*( V1X( 1 ) *V1 X( 1  ) 
))*PNiASS 
)  =  PKE1 
=  0.0 
=  O.C 
C   I  ■  1,6*4 

I  ♦  1 
C   J  =  IP  1,  17  2 

-1)   5.?r.  ,  5?r:  ,  51)0 


+  Vl  Yd  )  *VlY(  1  )  +  VIZ  (1  )  « 


-65)  540  ,53C,  550 

-  88)  60C  ,  6CC540 
R  1  X  (  J  )  -  Rl  X  (  I  ) 

=  A8SF(DRXJ  " 

RX  -  ROF  )  545,600,600 

R  1  Y  (  J  )  -  Rl Y (  I  ) 
=  APSF(CRY) 
RY  -  ROE  )  55*0, 6CC,  600 

Rl/  (J)  -  R1Z(I) 
=  ABSF(DRZ) 
RZ  -  ROE)  555,600,600 
=  (CRX*ORX  +  CRY* CRY  +  CRZ* CRZ) 
ST  -  R0E21  5  5£  ,.->00,6C0 
=  St.RTF(CIST) 

PCIF(CIST)  -  EPSTP 
)=PPE(  I)  +0.5*  POT 
)  =  PPE(J)  +  C.5  *  POT 

-  T.-'OT  +  POT 
KUE 

=  PPt (  1  )  -  0.5  »  POT  18 
PPE(C)  -C.5  *  P0T18 
C   1=2.64 

)  =  32.?16t268E-C8*(VlX(  I  )  «VlX(  I  )+Vl  YU  )*VlY(  I)  +V1  Z(  I  )*VlZ(  I  )  ) 
=  TPKE  +  PKP  (  I  ) 
NUE 

=  TPKE  +  PKE  1 
TPKl   +  TPOT 

=  APSF( (TE-EVJ/LV)  *1.CE2 

BLOCK  J 

WRITE  CO  1  PUT  I  APE  3,10 

WRITL  OUTPUT  TAPE  3,  12,  BX,  BZ,  EV,  CT I 

WRITE  OUTPUT  TAPE  3,  050 
S5C  FORMAT!  /,   19H  LATTICE  LOCATIONS  ,/) 

WRITE  OUTPUT  TAPE  3,860 
360  FORMA T(9CH    ATOM        CX         DY         CZ         VX  VY 

1  Ml  KE         PE       /) 

SHIFT  =  1.CE-C5 

CO  90 C  1=  1,6  4 

E  N I  ?  (  C  )  . 

CX  =  R  IX(  I  )-RXR(  I  ) 

+  SSK(CX),SLJ(L  +  2)  ,+SCM( 777778)  ,+FS8(SHlFT),AJP5(C65),  IM?(  1  )  . 

CY=  RlYl 1 )  -  RYR(  I  ) 

+  SSMCY)  ,SLJ(  L  +  2)  ,+SCM(77777B)  ,  + FSB  (  SH  I  F  T  )  ,  A  J  P  3  (  67  5  )  ,  IM2(  1  )  . 

CZ=  R  1Z(  I  )-  RZRU  ) 

+  SSK(CZ),SLJ(L+2)  ,+SCM< 777778) , +FSB( SHI  FT ) , A JP3 ( 895  )  , I N  I  2 (  1  )  . 
885  ENA2I  C)  ,  UPI900). 
890  WRITE  OUTPUT  TAPE  3 , 995 , I , CX , C Y , CZ  ,  V 1 X (  I  )  ,  V  1  Y (  I  )  ,  V 1 Z (  I  )  , P K E (  I  )  , 

1PPC(  1  ) 
90C  CONTINUE 
89  5   FOR  VAT { I10,3F10.5,3E1C.2,2F10.2) 


965 
6  75 
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HLCCK  K 


1CC5 


1C1C 


1015 


CCNV 
0T0C 
UTQRM 

i )  i  0  2  ft 

CTP?.," 
HT02R 
RX?  = 
RY?  = 

rso  =  ( 

RYl* 

XIP  =  R 
RXl  = 
VX2  = 
VY2  = 
VYl  = 
FX1  = 
FX2  = 
FY1  = 
FY?  = 
VXl  = 
INDEX 


C.C1 
CI  0C 

DT 

=  uf 

=  DTOR 

=  0 
0 


/DTI 

*  CPMV 

crw  «  ccnv 

£  i'"  «  tf!NV 
vi  *CONV 

T02RM1  *  CO\'V 


:-3.o 

R!  F( 


)MBZ-3.0)  +  (BX-3.C  )* (BX-3.0) 
DSQJ 


V 

1  = 
M 
C. 
C 
BZ 

sc 

Y1 

SCRTF{R0t2-l)SU  -  1.0E-06 

0.0 
0 
0 
0 

c. 
c. 
c. 


.c 

.c 
•  c 

0 

c 

c 

S.RTF(EV  *  C.C635U  /  P^SS)  *  5.S11811ECU 

c: 


CX  =  RXl-  RXl": 
CY  =  RYl-  RY2 
DIST  =  C/«DX  + 
IF (CI ST-R0E2) 


DY«OY 

101 C,  1010, 


BLOCK  LI 


1070 
BLOCK  L2 


CIST  =  Si,  RTF  (  01  ST  ) 

FORCE=  FfRCUMST)  -  LPST 

FXl  =  OX/  DIST  *  FCRCE 

FX2  =  -FXl 

FY1=   DY/DIST  *  FTRCt 

FY2  =  -  FY1 

IF(  INCEX)  101 5,  10  15,  10?0 


BLOCK  Ml 


RXU 

RYU 

PX2T 

RY2T 

FX1T 

FYl  T 

FX2T 

FY2T 

VX1T 

VYl  T 

VX2T 

VY2T  = 

VXl  * 

VY1  = 

vx:  = 

VY2  = 
RXl  = 
RYl  = 
RX2  = 
RY?  = 
INOEK 
GO  TO 


R  K  1 
RYl 
RX2 
RY/. 
FXl 
FYl 
■    FX2 
FY2 
Vxl 
VY  1 
'  VX2 
VY^ 
FXl 
F  YT 
F  X.. 
TY 

(VXl 
(VYl 
\V\2 
(  VY2 

1CC5 


i)TCRM  1 

DTORM 1 

OTORM 

D1GRM 

VX1T  ) 

VY1T  ) 

VX2T  ) 

VY2T  ) 


+   VXl 
+  VYl 
+  VX2 
+   VY2 
OTOC 

*  CTOC 

*  CTOC 
*uTOC 


■       RXU 
+       RY1T 
+    RX2T 
+    RY2T 
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BLOCK  V2 


1C3C 


VX  1  = 
VY1  = 
VX2  = 
VY2  = 
RX  1  = 
RY  1  - 
RX?  = 
RY2  = 
INDEX 
GU  TO 


(FX' 
{  F  Y  1 

(FX.: 

(FY. 
(  VaI 
(VYl 
(VX2 
(  V  Y  2 

Tcc5 


+  FX1T)«  0  1 C?RM 1  + 

+  FY1T  )*  CT02RM1 
+  FX2T) *  LTC2RM 
+  FY2T) *0TU2RM    + 
+  VX IT  )*  CTCC   + 


VY1T 
VX2T  ) 
VY2T.) 


DTOC 
«  DTOC 

*ltcc 


VX1T 

+  VYU 

VX2T 

VY2T 

RXlT 

RY1T 

RX2T 

RY2T 


107C 


10  75 

10800 
1 

65U 

655 

1 

660 

6650 
1 
2 

675 

3  COO 
1 

SIC 

1100 


P2KE1 
P2KE2 

TKE  = 

WRITE 

writ: 

WRl  TL 
WRITE 
FORMA 
WRI  TE 
FORMA 

F12. 
WRI  TE 
FCRMA 
WR  ITE 
FORMA 
PULLH 
WR  ITE 
FORMA 
WRI  TE 
FCRMA 

E  W/ 

LH 
WR  ITE 
FORMA 
WRITE 
FORMA 
PE 
WRITE 
FORMA 
WRITE 
FCRMA 
E  TWO 
INETI 


=  C.  : 

'  3c 
P2K 
OUT 
OCT 
CUI 
CUT 

T{// 
OUI 

T( 

5, 
OUT 

T(31 
CUT 

T(  ]u 

T 

our 

T  {  F  1 

OUT 
T(/, 

R  TO 
Y= 

OUT 
T  ( 

OUT 
T  ( 
R  ER 

OUT 
T(  3f 

OUT 
T(// 

PUO 
C  EN 


18  0U 

.916 

E  1  + 

PUT 

PUT 

PUT 

PUT 

//, 

PUT 

15 
UH 
PUT 
H  RE 
PUT 
H  TA 

KE 

PUT 
O.C, 
PUT 
1  UH 

LAT 
E9.2 
PUT 
/, 
PUT 
/, 
RCR 
PUT 
15.5 
PUT 
»  71 
Y  PR 
ERGY 


1C2E 

3263 
P2K 
TAPE 
TAPE 
TAPE 
TAPt 
1  CH 
TAPE 
H  Kl 
TO 

rAPE 

SCLT 

TAPt 

RGET 

T  A  R3 

I  APE 

1PE? 

TAPE 

PULL 

-TA} 

,  *H 

TAPE 

16H 

TAPE 


TAPE 
t  F 

TAPE 

H  CO 


PIE 
MAX 


-C8 

E-0 
E2 
3t 

-   » 

3, 

3, 
RES 

7 

N~E 
T  E 

3, 
S  0 

3, 

VA 
ET 

3, 
5.6 

3, 
ET 
GET 
I 

LAT 

7. 

oh' 


7 

ic! 
?, 

MPA 
M  , 
IMU 


BLOCK  N 

* ( VX1 *VX1  +  VYl  *V 
8  *(VX2*  VX2   + 


Yl ) *PMASS 
VY2»  VYt) 


15 
10 

12, 

107 
ULTS 
1080 
NRG 
NRG= 
t5M 
F  LA 

655 

LUES 

P 

o60. 

,  tie 

6  65, 
VALU 

=  F 
-   E 

675 
r  ICE 

000 

) 
810 
5  ) 

1  10 
RISC 
/,  6 


B  X  ,  B  Z  ,  E  V  ,  C  T  I 


OF  SIMPLE  TW 
t  P2KEltP2KL? 
PROJ=  ,  F12.5 

,  F12.5,  12H 


0  BOCY  INTERACTION    ,/  ) 
,TK2,EV 
,  18F     KIM  ENRG  TARG=  , 
ENRG  IN=  ,  F1C.1,  /////) 


TT ICE  INTERACTION  ,/,  ) 


E  STEPS    PCT  E  TARG 
/R  TC  LAl-^ULLET 
, PKE8R,  PEA3 


T  W/R  TO 


, //,9CH   TIM 
CT  £  TARGET  W 

CYCLE,  P0T1R 
.6,E25.6  ) 
MKE1,  C0L3,  V 
ES,//,15H    K 

7.1,25H     VELOCITY  IN  M/S EC   X=  1PES.2, 
9.2   ) 


1X(  1  )  ,Vl Y(l  ), VlZ (  1  ) 

IN  ENERGY  =  F  7.  1  , 


PCT 


VALUES    ) 

PCT  ENERGY 

,  TPOT  ,  TPKE 


KIN  ENERGY 
TE,  PERE 


TCT   ENERGY 


N  VALUES  OF  L 
5H  POTH  INTER 
) 


ATTICE  REACIICN  AGAINST  SIMPL 
ACTIONS  STOPPEC  WHEN  TARGtT  K 
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I  1C5 

II  1C 
1115 

1  1  16 

1  1  If 
1  1  19 


1  12C0 
1 
2 
3 

4 

11300 
1 

2 
3 
U 
5 

1,1  400 
1 


T2RA 

B2SA 

SRE8 

SRF1 

A  T  L  R  A 

IF(  AT 

TLRA 

CO  TO 

TLRA 

ARLSA 

IF(AB 

CLSA 

GC  TO 

PLSA 

RBMT 

ET2  = 

ETL  =  ( 

RET  =  E 

RKEP  = 

EABL 

RSA2L 

RRA2L 

ROT  1  = 

URITE 

FORMA 

18H 

/  KC 

KE  P 
S  = 
WRITE; 
FORMA 

39H 
TICE 

1  9F 
BODY 
LE  TW 
WRITE 
FORMA 
T  AT 
WRITE 


=  A 
=  A 
=  S 
=  _S 

LRA 

=  A 

1  1 

=  C 

LSA 

=  A 

1  1 

=  c 

=  p 

(P 

PKE 
TL/ 

PKE 
=  T 

=  B2 

PK 

CU 

T(/ 


I  ANF( 

TANF( 
..RTF( 
C  R  T  F  ( 
A  B  S  F  ( 
-1.0) 
COSF( 
15 
.C 

ABSF( 
-1  .0 
COSF< 
19 
.  0 

f<  ASS/ 
JKE2/ 
-R/EV 
LT2 
1/P2K 
c-PKE 
SA/BL 
1 2  R  A  / 
E  1/EV 
TPUT 
,  4  *H 


APSF  (VY2/VX2)  )     *57.296 

ABSF (VY1/VX1 ) )     *    57.296 

PKE8  R/32 . 91 63268 E-OB ) 

PKE1     /     (PMASS*C.518041C?E-08) ) 

VI  Y(  6  )  )/SREP 


1105,  1 
A  TL  R  A  ) 


10,1110 

*    57.296 


VI Y( 1 ) J/SRE1 
)     1116,1118,1116 
A"LSA)     *    57.296 


63.  54 
EV)*  1 .0E02 
)*1.CEC2 

Ll 

l-PKEfR-PPE(8)-PPE( 1 ) 
SA 

TLRA 


,SFER 
TWO 


TPA 
RfJ 

F  6  .  3 

OUTPUT 
T(/,33H 

SCATTE 
=   F7.2, 

/  SCATTE 
C  RCUY  / 

OUTPUT 
T(  /,  21 
END  OF  R 

OUTPUT 


TAPE 

P 

LATT 

TWH 
BODY 
) 

TAPE 

PO 

RING 

//,3 

LAT 
R  I  NG 

REC 
TAPE 
H 

UN    / 
TAPE 


3,  1 120,ET2,ETL,RET,RKEP,RRMT 
ERCENT    KL    TRAN'SPERED  .      TWO     BCCY=       F9.3, 

ICE  =    P9.3,//,     5cH  RATIO    KL     TRANSFcR     IN    LATTICE 

F.OCY    =       F7.5,//tU8H  RATIO    KE    PRCJ     IN    LATTICE     / 

=    F11.5,//,     36H  RATIO    GULLET    MASS    /     TARGET    MAS 


5, 
T     + 

AN 
5K 
TIC 

AN 
CIL 

7 


1    1 

K 
GL 


GL 
A 

,1  1 

IMPAC 

ORIG 

3,15 


30,EABL,p2SA,3LSA, T2RA, TLRA, RSA:l,°RA2L 

E     ABSORBED    FY    LATTICE    =       F12.3     ,//, 

E  TWO  BODY  =  F7.2,  1 9H  LAT 

RECOIL  ANGLE  TWO  BODY  =  F7.2, 

=   F7.2,//,63H    RATIO  SCATTERING  ANCLE  TWC 


E  LATTICE  =    F9.5,//,55H 
NGLE  LATTICE  =   F9.5     ) 
1*0,  XI  P,  RBFI 
T  PARAMETER  =   F7.4,   S4H 
INAL  ENERGY  =   1PE13.5 


RATIO  RECOIL  ANG 


RATIO  KE  BLLLE 
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BLOCK  R 


116C 

1  17C 

1  18C 


C  (  1  ,  L)  =  f  X 
C ( 2 , L  )  =  HZ 
C{ 5,L  )  =  X  IP 

IF ( P2-3.0)  1  160,1  ItC,  1 17C 
C(U,  L)  =  C.C 
1180 
=  ACCSF  (  (2 


go  tu 

C(U,L) 
C  (  c.  ,  L  ) 
C(6,L  )  = 
C(7,L  )  - 
C(8,L  )  = 
C(9,L  )  = 
C( 10, L)  = 
C(  1  1,L 
C(  1?, 
C(  1  3, 
C(  1U, 
C(  15, 
C(  U, 
C(  17, 
C(  IP, 
C  (  1  9  , 
C(20, 
C(21, 
L  =  L+1 


l-l 
tT2* 1 .OE-02 

P2SA 

ETL»  1  .0E-C2 

BLSA 

RET 

- PPE(8 )  +  PKE8R 
=  TLRA 
L)*PKE  1 
L)=P2KE1 

PKE8R 

P2KE2 

PERE 

YCLL 

PPEt  1  ) 

PPF_(  8  ) 

I"2RA 


0)/XIP)»  57.296 


L)  = 
L)  = 

)  = 

)=C 

)  = 

)  = 

)  = 


)=RlY( 1 ) 


BLOCK  T 


r- 

,  3  , 


BZ  =  BZ-C. 1 
H  (BZS-B2  >5 

1?0C  CONTINUE 

IF(SENSE  SWITCH 

UOOC  CONTINUE 


1  ZOO 
1  ) 


UOCC 
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BLOCK  W 


13  5C 


13  60 
1400 


1*4  50 
150C 


1510 
1530 


LMAX 
WR  I  T 
WRIT 
FCRV 
1SANG 

00  T 

wr  i  r 

FOR^ 
DC  1 
WRIT 
J  = 

1  = 
WRIT 
WRIT 
WRIT 
WR  I  T 
FORM 

1  ICE 
2ICC 

WR  I 
1  (C( 
2(C( 
3(C( 

DO 

WRI 

J=l 

1  =  1 

WRI 

WRI 

WRI 

WRI 

FORM 
1  T  ICE 
2F00Y 


=  L- 
L  GDI 

c  rui 

AT(  // 


?l. 


0  J 
E  GUI 
AT(3F 

SCO 
E  0U1 
12 
14 

E  OUT 
E  HUT 

E  cur 
l  our 

AT(  // 
//,3. 
) 

e  cur 

,L), 
»Ui 
,L), 
53C 

E  CU1 


1 

PUT  TAPE 
PUT  TAPE 
,  1  1  0  H 

TW(L) 
■  1  ,  LMAX 
PUT  TAPE 
9.3,  no. 
M  =  1  » 5 
PUT  TAPE 


3,10 

3,1350 
x 
S ANGLE 

5,  U00, 
2,F1C4, 


Z        B 
TM(L)/TM(2) 


ANGLE     TM2) 
TE  TAR     R ANGLE    i/  ) 


(C(  I,  J  ),  1=  1,11   ) 

F10.2,F10.4,Fl0.2,F10.4,F13.4,FlC.2) 


3,15 


PUT  TAPE 
PUT  TAPL 
PUT  TAPE 
PUT  TAPE 
, 33X53H 
X57H    ROT 

PUT  TAPE 
L=12,20) 
L  =  21  ,27) 
L=33,35) 
M  =1,5 
PUT    TAPE 


5,  1U2  5 

3,  10 

3,  1430,EV,DTI 

3,14  50 
TCP    NO.     IS 
TOW    NO.     IS 


KINETIC 
KINETIC 


ENERGY 
ENERGY 


IN 
Ii\ 


EV 
EV 


CF 
OF 


BULLET 
TARGET 


IN 
IN 


LA  TT 
LATT 


,1050, 
(  C  (  I  ,  L 
(C(  J,L 
(C (  I  ,  L 


(C(J,L),  L=1 
)  ,  1  =  12,20), 
),  L=28,32), 
),  L=33,35), 


,11),  (C(  I,L)  ,  L=l,  1  1  ), 
(C(J,L),  L=2  1,J.  7), 
(C(I,L),  L =2 0,32), 
C( J,  36)  ,C(I,36) 


3,1 


1540 
ft  50 


WRI 

1  (C( 

2(C( 

3(C( 
DC 
WRI 
J  = 
I  = 
WRIT 
WRIT 
WRI  T 
WRI  T 
FORP 

1  T  I  C  £ 

2C  BC 


F  CUT 

L  CUT 
E  CUT 
E  CUT 
A  T  (  /  / 

,//,: 

PROH 
E  CUT 
•  L), 
,L)  , 
,1), 
55C  P 
E  CUT 
1  4 
1  5 


WRIT 

1  (C(  J 

2(C(  I 

3{C(  J 

WR  IT 

J=17 

1  =  16 

WRIT 

WRIT 

WR  IT 

WRIT 

15  60  FORM 

1  NO. 

WR  I  T 

1  (C(  J 

2(C(  I 

3(C(  J 


E  CUT 

e  cur 

E  CUT 
E  CLT 
ATI  // 
,//, 
UY  PR 
E  CUl 
,L)  , 
,1)  , 
,L), 

e  cur 


E  OLl 
E  CUT 
E  CU1 
E  0U1 
M(// 
IS  P 
E  OUT 
,L), 
,L), 
,L), 


PUT 
PUT 
PUT 
PUT 
,  33 
0X64 
LEM 
PUT 
L=  1  2 
L  =  21 
L  =  3  3 
=  1, 
PUT 


PUT 

PUT 

PUT 

PUT 

,  33 

27X6 

OBLC 

PUT 

L=12 

L  =  21 

L  =  33 

PUT 


TAPE 
TAPE 
TAPE 
TAPE 
X55H 
H    BO 

) 
TAPF 
,20) 
,27) 
,35) 
5 
TAPE 


1  APE 
TAPE 
TAPE 
TAPE 
X53H 
6H  B 
M  ) 
TAPE 
,2  0) 
,27) 
,35) 
TAPE 


,14: 

10 


3,  1430, EV,CTI 


IS 
IS 


.-.,1510 
TCP  NO. 
TTCM  NO. 

3,  1850, 
,  (C( I ,L 
,  (C( J,L 
,  (C( I,L 

3,15 


,1425 

,1C 

, 1U30,EV,DTI 


KINET  IC 
KINETIC 


ENERGY 
ENERGY 


IN 
IN 


EV 
EV 


CF 
CF 


BULLET 
BULLET 


IN 

IN 


LAT 
TWO 


(C(J,L),  L=l 
),  L=12,20), 
),  L=28,32), 
),  L=33,35), 


,11),  (C( I, L)  ,  L  =  l  ,1  1  ) 
(C( J,L) ,  L=2 1 ,27) , 
(C(I,L),  L=2  8o2), 
C( J, 36)  ,C(  1,36) 


-, 1 540 
TCP  NO. 
CTTOM  NO 

3,1850, 


IS  KINLTIC  ENERGY  IN  EV  CF  TARGET  IN  LAT 
IS  KINETIC  ENERGY  IN  EV  OF  TARGET  IN  TW 


(C(  J,L 
(C(  I,L 
,15 


(C( J,L)  ,  L=1  , 
),  L=12,20), 
),  L=28,32), 
),  L=33,3S), 


L=l,  1  1  ) 


1 1 ) ,  (C(I,L) ,  L=l, 
(C(J,L),  L  =2  1,27), 
(C(  I  ,L)  ,  L=28o2  ), 
C( J, 36)  ,C(  I, 36) 


PUT 
PUT 
PUT 
PUT 
,  4 
ERC 
PUT 
L=l 
L  =  2 
L  =  3 


TAPE 

TAPE 

TAPE 

TARE 

4X32H 

ENT    E 

TAPE 

2,20) 

1  ,27) 

3,35) 


3,1425 
3,10 
5, 1430, 
3, 1560 
TCP  NO. 
RRCR   ) 

3,1850, 
,  <C(  I,L 
,  (C(J,L 
,  (C(  I,L 


EV,DTI 
IS  NUMBER  OF  TIME  S TEPS, / / , 46X28F  BOTTCM 


(C( J,L) ,  L=l 
),  L=12,20), 
),  L»28,32), 
),  L=33,35), 


,11),  (C( I ,L) ,  L=l, 1 1 ), 
(C( J,L)  ,  L  =  2  1,27), 
(C(I,L),  L*28.32li 
C( J, 36)  ,C( I,  36) 
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157C 
15RC 


1590 

1600 


1610 
1620 


{1630 

•1640 


1425 
1430 

1350 


un  is 

WRITE 
J  =  R 
I=«s 

WRI  Tl 
WRITE 
WRITE 
WRI  TE 
FORMA 

1  BOTT 
WRITE 

1 (C( J, 

2  (  C  (  I , 
3(C( J, 

CO  16 
WRI  Tb" 
J=l  1 
1=20 
writ: 
wri  re 

WRITE 
WRITE 
FORMA 

l  roiV  n 

WRITE 

1  (C( J, 
2(C(  I, 
3(C( J, 

DO  16 
WRI  Tf- 
J=19 
1=18 
WRITE 
WRI  TE 
WRITE 
WRITE 
FORMA 
1 Tl  ICE 

2  LATT 
WRITL 

1 (C( J, 

2(C( 

3(C( 
CO 
WRI 
J  =  2  1 
1=12 
WRI  TE 
WRITE 
WRITE 
WRITE 
FORMA 

1 ,32X5 
WRITE 

1 (C( J, 

2(C{  I, 

3(C( J, 
FORMA 
FORMA 

1EP  MU 
FORMA 

1  « 

2 

3 


P.C       M=  1,2 
OUTPUT  TAPE  3,15 


CUTPU 

cur  pu 

CLTPU 
OUTPU 

CM  '.0. 
PUTPU 
L),  L  = 
L),  L  = 
L),  L  = 
CO  h  = 

oui  pu 


t  tap: 

T  TAPE 
T  TAPE 
T  TAPE 
UOXi'JH 
IS  sc 
T  TAPE 
,20) 
,27) 
25) 


12 
21 
3  3 
It* 

T  TAPE  3,15 


M 
3,1 

3,  1 

2ii 

TCP 

att: 

5,1 
,  (C 
,  (C 
,  (C 


425 

c 

4*0, LV, DTI 
570 
NO. 

RING 

BSO, 
(  I,L 
(J,L 
(  I,L 


IS  SCATTERING  ANCLE  IN  L A T I  I CE , / / , 34 x52H 
ANGLE  OF  TWO  BOCY  PROBLEM 


(C( 
), 
), 
)t 


J,L),  L=l 
L=12,20) , 
L=23,32), 


L=33,35) 


) 

11),  ( C ( I , L )  ,  L  = 1 
(C( J,L) ,  L=2  1 ,2  7) 
(C(  I,L)  ,  L=28,  -i2) 
C( J, 36)  ,C(I, 36) 


1  1  ) 


OU 
CU 
OU 
OU 
T(/ 
C. 

OU 
L)  , 
L), 
L), 
20 
OU 


TPUT 

1  PUT 

TPUT 

TPUT 

/,  4 

IS  R 

TPUT 

L=l 

L  =  2 

L  =  3 

M=1 

IPUT 


tap: 

TAPE 
TAPE 
TAPE 

:*35h 

tCOIL 
TAPE 
2,20) 
1  ,27) 
S35) 
•  5 

TAPE 


3,14 

3,10 

3,  U 

3,15 

TOP 

ANGL 

3,18 

(C( 

(C( 

(C( 

3,15 


30,tV,DTI 

00 

NO. 

E  r 

50, 

I.L 

J,L 

IiL 


LATTICE,//, 36X47H  ROT 

) 

(I  ,L)  ,  L=l,  1  1  ), 


,  IS  RECOI.L  ANGLE  IN  LATTICE,//, 
)F  TWO  BODY  PROBLEM   ) 

(C( J,L) ,  L  =  l ,1  1  ),  (C(I  ,L)  .  L=l, 
.),  L=12,2C),  (C(J,L),  L=21,27), 
.),    L=20,32),     (C(1,L),     L=28,32), 

),    L=33,35),    C( J, 36)  ,C( I, 36) 


I, 
J, 
16 

TE 


OUT 
OUI 
OUT 
OUT 

T(// 

,//, 

ICE 
OUT 

L), 

L), 

L), 

140      W 

OUT 


GUI 
OUT 
OUT 
OUT 
T(// 
6H  B 
OUT 
L), 
), 
), 
(// 
(/, 
LTIP 
T(// 


«  t/ 
3X9F10.3 


PUT 
PUT 
PUT 
PUT 
,32X 
31X 
) 
PUT 
L=12 
L  =  21 
L=33 
=  1,2 
PUT 


PUT 

PUT 

PUT 

PUT 

,  36X 

OTTO 

PUT 

L=12 

L  =  21 

L=33 

//// 

ex2 

LIcR 
//,3 

« 
,3X1 

» 

,/// 

* 


3,  1U2  5 
3,10 

3,  1430,  LV, DTI 
3,1610 
TCP  NO.  IS  POTENTIAL  ENERGY  IN  EV  OF  TARGET  IN  LA 


TAPE 
,20) 


7) 
35) 


3,1850, (C( J, L)  ,  L=1 
(C(  I,L)  ,  L  =  12,20) , 
(C( J,L)  ,  L=28,32) , 
(C(I,L),  L=33,35), 


,11),  (C(I  ,L)  ,  L=l 
(C( J,L)  ,  L=2  1,27) 
(C(I,L),  L=28,32) 
C( J, 36)  ,C(I, 36) 


1  1  ) 


TAPE 


T  A  PE 

TAPt 

TAPE 

TAPE 

U7H 

M  N3 

TAPE 

,20) 

,27) 

,35) 

///) 

5H  R 


-' , 


15 


3,1 
3,1 
TCP 
.  IS 
3,1 
,  tC 
,  (C 
t  (C 


425 

C 

4*0, FV, DTI 

6  50 

NO. 


KIN 
8  50, 
(I,L 
(J,L 
(I,L 


IS  PULLET  PENETRATION  IN  LATTICE  UNITS,// 


ETIC  ENERGY 
(C(J,L),  L=l 
),  L=12,2C), 
),  L=28,32), 
),  L=33,35), 


IN  tV  OF  BULLET  IN 
,11),  (C(I,L)  .  L=l  , 
(C( J,L)  ,  L=2  1  ,27), 
(C(I,L),  L=2  8,.2), 
C( J, 36) ,C(I, 36) 


LATT  ICE 
1  1  ) 


) 


XI  IF 


LLLET  INCOMING  ENERGY  =  ,F8.1,3lH   t_V      TIME  ST 
F5.3,  34H      BULLET  INCIDENT  CN  (100)  FACE   ) 
1  C  .  3  ,  /  ,  1  1  1  H  «  »  * 

1F10.  3, //////, 1 3X9F10.3, /,  18X82F  *  * 

ik-  *  *  -ft  *    / 

///,  23X7F10.3,/,28X62H    *  *  * 

«  *,/, ? 3X7F1C. 3, //////, 33X5F1C. 3,  / ticX 

*  *  »  *,/, 33X5F 10.3,//////  ,U 

8X211    *,/,53XFlC3) 


73X3F1C.3,/,»*8X22! 
B3,/,5 

END 


*,/,U3x3F10. 3, //////, 53XF 10. 


88 


BLOCK  Y 

FUNCTION  POTF(OIST) 

COMMON  E1FP.AFAC, ETFF 

POTF  =  ETFP«(l.O-Kr.)IST*AFAC)»»0.803U)«*{-3.7BU)/DIST 

RETURN 

END 

BLOCK  I 


FUNCT 

COMMO 

Y=(CI 

YP=1. 

Z  =  YP* 

FORC  = 

RETUR 

END 

END 


ICN  FORCU,  1ST) 
IS  El FP,AFAC, ETFF 
ST«4FAC)  **C.  8C3U 
C  +  Y 

ETFF*Z»( 1 .  0*2.  9999*  Y/YP)/(DIST»DI  ST) 
N 
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APPENDIX  VI 
POTENTIALS 
1.   Born-Mayer 
A.   Potential 


4>M  »  Ae 


r  ■  nearest  neighbor  separation  at  zero  pressure  and  absolute  zero 
r  for  copper  =  2.51A 

The  three  different  sets  of  constants  used  in  the  potential  form  are 
identical  to  those  used  by  Gibson  et.  al. 

Potential  A(ev)  /O 


1  0.0392  16.97 

2  0.051  13.00 

3  0.1004  10.34 

Potential  No.  2  was  used  extensively  in  the  program.   It  can  be  reduced 

to  a  simplified  form  in  the  following  manner: 

-  13.©  (  r- ass  I  W  55, 
9(r)s  O.OSI  a  e 

if  distances  are  expressed  in  lattice  units 

(r)m    0.OSI  *  e 


4> 


in  the  potential  function  subroutine  this  becomes 

P0TF  =  0.051  *  EXPF(13.0-DIST  *  9.208) 
Potentials  one  and  three  become; 

1)  POTF  =  0.0392  *  EXPF(16.97-DIST  *  8.9996) 
3)   POTF  =  0.1004  *  EXPF(10.34-DIST  *  7.3115) 
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B.   Force 


=  -Ae 


-P(r-r.)/ 


/° 


substitution  of  constants  for  Potential  No.  2 

Tr         P.  55 1 

if  distances  are  expressed  in  lattice  units  and  electron  volts  changed  to 
joules; 

^rdfrfl  .   "*  0*.o-S.2o8r) 

<*  L^r)J     .    -o.4iW4*io     x    e 

d  r 

in  the  force  function  subroutine  this  becomes; 

FORC  =  0.41614E-09  *  EXPF(13.0-9. 208*DIST) 
forces  for  potentials  one  and  three  become; 

1)  FORC  =  0.3127E-09*EXPF(16.97-8.9996*DIST) 

2)  FORC  =  0.65C7E-09*EXPF(10.34-7.3115*DIST) 


2.   Thomas-Fermi-Firsov 
A.  Potential 


V(r)     =  1   E 


TF 


_IE 

r  , 


*(£)=+w 


Thomas-Fermi  screening  function 


-<*: 


4W- 


i  +  <£, 


fc 
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•d*  -  o.  8034 
<**  5  3.731 


■TF 


"  **'*'*%„ 


i 


QTP  = 


(45T  +  ^  f 3 


first  Bohr  radius  of  hydrogen 


K*l 


-ll 


£^  =  21^233^^      qTp=.qioO(oa^o    m 


*TF 


let  x  ■  r/a 


TF 


and  ATFR  =  1/a    (in  lattice  units) 
TF 


then  x  =  ATFR  *  DIST    (DIST  in  lattice  units) 

♦«•[>*  (A  rr  *  H 


tfzl 


_p6 


let  ATFR^       =  AFAC 


4>(X)  =    h  +  C^F/qc-Disr) 


0,SO3H  -i 


-3.13^ 


V(r)  = 


•TF 


J.  $(x)  =    Etp  ;  foOO 

X  2  *  ATFR-  DIST 


let      ETFP   =  E      /    2  *  ATFR 
then 


POTF  s     ETFP 

DJST" 


o.eo&i : 
\.o  +  (7}FAC*D1ST) 


3^34 
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B.      Force 


set     x  =  r/a 


TF 


^  (VQ)) 


^m) 


s* 


TF 


^r  2 


y2  X        dx 


'TF 


set 


Y-(%.) 


then 


di  ( <H*>) 


3      X     UY 


^(vlCo)    r    Etc  «Qtt 


4>W 


if  the  following  substitutions  and  appropriate  changes  in  units  are  made 


where 


pIST  *  "DIST 


ETFF  =  ETF  *  aTF  /2 


Y^ 


x  =  DIST  *  AFAC  (DIST  in  lattice  units) 

Y=  (^)oS034  YP  =  »-o  +  Y 

* = (YP) 

ETFF  and  AFAC  are  defined  at  the  beginning  of  the  program  -  Y,  YP  and  Z 
are  calculated  in  the  subroutine  each  time  it  is  used. 
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